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FOREWORD 


The  research  results  contained  in  this  technical  report  were  performed  under  an  Air 
Force  Contract  F33615-91-C-3205  entitled,  “Finite  Element  Thermal  Acoustic  Nonlinear 
Analysis  of  Aerospace  Structures.”  The  report  is  based  essentially  on  a  Ph.  D.  dissertation 
prepared  by  Roger  R.  Chen  under  the  supervision  of  Dr.  Chuh  Mei.  The  document 
presents  a  finite  element  formulation  and  solution  procedure  for  the  prediction  of  nonlinear 
response  of  composite  plates  of  arbitrary  shape  subjected  to  intense  acoustic  and  thermal 
loads.  The  work  was  conducted  at  the  Department  of  Aerospace  Engineering,  Old 
Dominion  University.  Mr.  Kenneth  R.  Wentz  and  Dr.  Howard  F.  Wolfe,  WL/FIBGD, 
Structural  Dynamics  Branch,  Flight  Dynamics  Directorate,  were  the  technical  monitors. 


Chapter  1 


INTRODUCTION 

1.1  Background 

“Jet  Crash  Off  Italy  Kills  35”  The  New  York  Times  reported  on  January  11,  1954. 
“Rome,  Jan.  10  —  Thirty-five  persons  were  almost  certainly  killed  when  a  British  Comet 
jet  airliner  crashed  into  the  sea  this  morning  about  halfway  between  the  islands  of  Elba 
and  Monte  Cristo,  off  the  Italian  western  coast.”  It  was  Monday,  when  people  went  to 
work,  they  read  this  news  very  sadly.  The  jet  airliner  was  the  Comet,  the  first  propelled 
by  jet  turbo  engine. 

The  news  surprised  the  world  aircraft  designers  and  manufacturers.  After  exhaustive 
investigation  and  tests  on  components  of  the  Comet,  it  was  concluded  that  the  accident 
was  caused  by  fatigue  failure  of  the  pressurized  cabin.  The  small  fatigue  crack  origioated 
from  a  comer  of  an  opening  in  file  fuselage.  Since  then,  extra  additional  attention  has 
been  focused  on  airframe  fatigue  design.  It  has  been  realized  fiiat  such  failures  can 
substantially  increase  the  maintenance  burden  and  life  cycle  cost  of  the  aircraft. 

In  the  late  1950’s,  incidents  were  reported  in  which  stmctural  components  close 
to  high  intensity  jet  exhausts  were  reaching  such  high  levels  of  vibration  response, 
due  to  acoustic  excitation,  that  fatigue  cracks  could  develop  and  spread  quite  rapidly. 
These  incidents  alerted  industry,  university  and  research  centers  to  the  possibility  of 
serious  design  problems  as  the  performance  of  aircraft  and  engines  increased.  Since 
then,  the  acoustically  induced  fatigue  failures  in  aircraft  have  been  one  of  the  major 
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design  considerations.  A  considerable  number  of  investigations  sponsored  by  USAF 
and  AGARD  have  been  carried  out  during  the  1960’s  and  early  1970’s.  This  series  of 
research  lead  to  design  monographs  of  every  type  of  metal  structures.  It  made  it  possible 
for  aircraft  designers  to  produce  satisfactory  structures. 

This  was  the  situation  towards  the  middle  of  1970’s.  Although  the  power  of  the  jet 
engines  was  still  increasing  dramatically,  the  pressure  levels  were  not  increasing  because 
of  the  use  of  the  higher  bypass  configurations  in  engine  design  which  are  needed  to  reduce 
the  community  noise.  This  lead  to  a  great  reduction  in  the  research  and  development 
activity  on  sonic  fatigue. 

In  the  mid  and  late  1980s,  new  interest  in  random  vibration  was  raised  primarily 
as  a  result  of  advances  in  high  speed  flight  The  missions  for  these  high  speed  flight 
vehicles  will  expose  structures  to  severe  acoustic  and  thermal  loads.  These  include  the 
F-22  advanced  tactical  fighter  (ATF),  the  supersonic  advanced  short  take-off  and  vertical 
landing  aircraft  (ASTOVL),  the  national  aerospace  plane  (NASP)  and  the  high  speed  civil 
transport  (HSCT).  The  new  vectored  thrust  propulsion  systems  on  ATF  and  ASTOVL 
provide  short  take-off  and  landing  capability,  and  also  increased  maneuverability.  Most 
designs  feature  engine  exhaust  locations  that  are  positioned  near  the  aircraft  mass  center, 
and  exhaust  jets  are  directed  either  onto,  or  nearly  onto,  the  aircraft  structure.  Estimates 
of  acoustic  loads  indicate  that  for  some  vectored  thrust  directions  most  of  the  aircraft  is 
immersed  in  an  acoustic  field  with  levels  well  above  150  dB,  with  levels  much  higher 
near  the  nozzle  (Mixson,  1988).  Exhaust  temperatures  may  exceed  1000°F  in  the  region 
of  nozzle,  and  therefore  the  structure  must  withstand  high  thermal  loads.  In  addition, 
the  HSCT,  ATF  and  NASP  will  fly  at  supersonic/hypersonic  speeds  and  will  be  exposed 
to  intense  in-flight  acoustic  and  thermal  environments  (Pozefsky  et  al.,  1989).  Due  to 
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aerodynamic  heating,  the  structures  will  experience  high  temperatures  with  large  thermal 
gradients. 

To  meet  increased  performance  requirements,  new  complex,  lightweight  structures 
and  advanced  materials  will  be  required.  The  complex  structures  under  consideration 
have  significant  uncertainties  in  fatigue  behavior  due  to  intense  acoustic  loads  in  the 
presence  of  high  temperatures.  The  thermal  environment  can  affect  acoustic  fatigue  by 
introducing  thermal  inplane  forces  and  thermal  bending  moments,  as  well  as  altering 
(temperature  dependent)  material  properties.  Such  thermal  effects  may  also  introduce 
large  distortions  and  snap-flnough  (or  oil-canning)  behavior,  alter  buckling  loads  and 
modify  vibration  characteristics.  The  intense  acoustic  loads  can  affect  fatigue  life  by 
introducing  large  deflection  geometrical  nonlinearities,  modal  coupling  and  multiple¬ 
mode  participation  (  Mei  and  Wolfe,  1986).  Such  high  sound  pressure  levels  may  even 
drive  the  structures  to  have  damping  nonlinearity  (  Mei  and  Prasad,  1987).  Because  of 
high  costs  and  difficulties  with  instrumentation  in  experiments  at  high  acoustic  intensity 
and  elevated  temperatures,  reliable  experimental  data  is  difficult  to  acquire.  Thus,  in  the 
design  process,  greater  emphasis  will  be  placed  on  analytical  and  computational  methods. 
This  brings  a  tremendous  challenge  to  the  analysts  for  predicting  nonlinear  response  of 
complex  structures  subjected  to  acoustic  and  thermal  loads. 

A  fundamental  challenge  of  thermo-acoustic  random  response  of  aerospace  structures 
is  the  multi-disciplinary  nature  of  the  problem.  The  thermal  environment  strongly  affects 
the  structural  random  response  because  of  the  development  of  restraint  forces  due  to 
thermal  expansion  and  the  change  of  material  properties  at  elevated  temperatures.  Since 
the  structural  response  due  to  high  levels  of  acoustic  loads  is  highly  nonlinear  and  strongly 
dependent  on  the  thermal  effects,  the  problems  are  thus  inherently  multi-disciplinary  and 
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nonlinear.  It  is,  therefore,  the  purpose  of  this  study  to  develop  an  analytical  formulation  to 
determine  the  nonlinear  random  response  of  composite  laminates  to  combined  acoustic 
and  thermal  loads  applied  simultaneously. 

1.2  Literature  Survey 

A  limited  number  of  papers  studied  the  thermal  postbuckling  of  laminated  composite 
plates.  Noor  and  Peters  (1983)  have  investigated  the  bifurcation  buckling  and  postbuck¬ 
ling  response  of  composite  plates  subjected  to  combined  axial  compression  and  uniform 
temperature  distribution  by  the  multiple-parameter  reduced  based  technique.  Recently, 
they  (Noor  et  al.,  1992)  have  studied  the  thermo-mechanical  buckling  and  postbuckling 
response  of  composite  plates  subjected  to  combined  axial  and  thermal  loadings.  The 
analysis  is  based  on  a  first  order  and  third-order  shear  deformation,  von  Karman  type 
of  nonlinear  plate  theory.  A  mixed  formulation  is  used  with  the  fundamental  unknowns 
consisting  of  the  generalized  displacements  and  the  stress  resultants  of  the  plate.  An 
efficient  multiple-parameter  reduction  method  is  used  in  conjunction  with  mixed  finite 
element  models.  Sensitivity  derivatives  are  evaluated  and  used  to  study  the  sensitivity 
of  the  postbuckling  response  to  variations  in  the  different  lamination  and  material  pa¬ 
rameters  of  the  plate.  In  a  paper  presented  at  the  35th  SDM  conference  (Noor  et  al., 
1994),  a  similar  study  of  the  composite  plate  with  cutout  was  conducted.  Huang  and 
Tauchert  (1986)  used  analytical  continuum  approach  and  studied  the  thermal  buckling 
and  postbuckling  behavior  of  simply-supported  antisymmetric  angle-ply  plates  subjected 
to  uniform  temperature  change.  Their  results  illustrated  the  effects  of  the  number  of 
layers,  the  ply  angles  and  the  aspect  ratio  of  the  plate  upon  the  thermoelastic  response. 
Chen  and  Chen  (1989,  1991)  have  studied  the  thermal  postbuckling  behaviors  of  lami¬ 
nated  rectangular,  antisymmetric  angle-ply  composite  plates  subjected  to  a  nonuniform 
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temperature  field  by  the  finite  element  mefiiod.  Based  on  the  principle  of  minimum  po¬ 
tential  energy,  the  nonlinear  stiffness  matrix  and  geometric  stiffness  matrix  are  derived. 
Their  results  reveal  that  the  thermal  postbuckling  behavior  of  composite  laminated  plates 
is  influenced  by  lamination  angles,  plate  aspect  ratio,  modulus  ratio  and  the  number  of 
layers.  Their  results  also  revealed  that  the  effect  of  temperature-dependent  mechanical 
properties  on  the  thermal  postbuckling  behavior  is  significant.  Librescu  et  al.  (Librescu 
and  Souza,  1991  and  Librescu  et  al.  1994)  recently  studied  the  static  postbuckling  of 
simply  supported  flat  panels  exposed  to  a  stationary  nonuniform  temperature  field  and 
subjected  to  a  system  of  subcritical  in-plane  compressive  edge  loads.  The  study  is  per¬ 
formed  within  a  refined  theory  of  composite  laminated  plates  incorporating  the  effect  of 
transverse  shear  and  the  geometric  nonlinearities.  Meyers  and  Hyer  (1992)  have  analyti¬ 
cally  studied  thermal  buckling  and  postbuckling  of  simply  supported  symmetric  composite 
laminates  under  uniform  temperature  change  using  flie  Rayleigh-Ritz  method.  Birman 
and  Bert  (1993)  investigated  the  effects  of  temperature  on  buckling  and  postbuckling 
behavior  of  reinforced  and  unstiffened  composite  plates  and  cylindrical  shells.  First,  the 
equilibrium  equations  are  formulated  for  a  shell  subjected  to  the  simultaneous  action  of  a 
thermal  field  and  an  axial  loading.  These  equations  are  used  to  predict  a  general  form  of 
the  algebraic  equations  describing  the  postbuckling  response  of  a  shell.  Conditions  for 
the  snap-through  of  a  shell  subjected  to  thermo-mechanical  loading  are  formulated.  The 
theory  was  also  applied  to  prediction  the  postbuckling  response  of  flat  large-aspect-ratio 
panels  reinforced  in  the  direction  of  their  short  edges. 

There  are  many  survey  articles  which  review  analysis  techniques  and  experiments  of 
nonlinear  random  structural  responses.  Mei  and  Wolfe  (1986)  have  presented  a  discussion 
on  analytical  and  experimental  techniques  to  predict  the  acoustic  fatigue  life  of  aircraft 
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structures.  They  discussed  die  problem  and  the  steps  taken  to  solve  the  problem  and 
reviewed  the  analytical  approaches  to  single  degree-of-freedom  (SDOF)  and  multi-degree- 
of-freedom  (MDOF)  linear  and  nonlinear  systems  under  random  excitations.  They  also 
reported  the  advances  of  analytical  prediction  and  experiments.  To  (1987)  has  also 
presented  a  comprehensive  survey  paper  on  the  analysis  of  nonlinear  systems  subjected 
to  random  excitation.  Methods  reported  to  be  applicable  to  both  SDOF  and  MDOF 
systems  include  equivalent  linearization  techniques  (EL),  the  Fokker-Planck-Kolmogorov 
equation  (FPK)  and  moment  approaches.  Clarkson  (1994)  has  recendy  presented  a  very 
comprehensive  sonic  fatigue  technology  review  report.  He  reported  that;  “  From  the 
early- 1960s  until  the  mid-1980s,  there  was  very  little  theoretical  development  for  sonic 
fatigue  prediction.”  The  design  monographs  for  most  common  aircraft  structures  were 
made  based  on  simple  theoretical  models  and  results  of  specially  designed  tests.  The  use 
of  advanced  composites  in  the  1980s  generated  an  increasing  interest  in  development  of 
more  sophisticated  theoretical  models,  because  the  much  wider  range  of  parameters  of 
composite  panels  made  creation  of  nomograph  based  on  tests  not  possible. 

A  limited  amount  of  investigations  on  structural  response  subjected  to  intense  acoustic 
and  thermal  loads  exists  in  the  literature.  Seide  and  Adami  (1983)  were  the  first  who 
studied  large  deflection  random  response  of  a  thermally  buckled  simply  supported  beam. 
The  thermal  load  and  the  acoustic  pressure  are  thus  considered  to  be  applied  in  sequence. 
The  well-known  classic  Woinowsky-Krieger  large  amplitude  beam  vibration  equation  is 
used.  The  Galerkin’s  method  and  time  domain  numerical  simulation  are  then  applied  to 
obtain  random  response. 

The  papers  by  Mei  and  Prasad  (1987  and  1989)  aim  to  explain  the  observed 
broadening  of  the  response  peak  and  its  increase  in  frequency  by  including  nonlinear 
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damping  as  well  as  large  amplitude  displacements  in  the  theory.  This  is  a  very  valuable 
formulation  because  damping  in  inherently  nonlinear  and  its  behavior  and  magnitude  is 
one  of  the  major  unknowns  in  the  work  to  date.  In  their  work  a  single  mode  analysis  is 
used  and  the  results  show  the  expected  broadening  and  increase  in  frequency. 

Most  recently,  the  Galerkin/numerical  simulation  approach  was  applied  to  simply 
supported  metal  and  orthotropic  composite  rectangular  plates  by  Vaicaitis  and  Arnold 
(1990)  and  Vaicaitis  (1991),  the  thermal  and  acoustic  loads  are  considered  to  be  applied 
simultaneously.  The  classic  von  Karman  large  deflection  plate  equations  including 
temperature  and  orthotropic  property  effects  are  employed.  The  thermal  effects  on 
rectangular  isotropic  plate  random  response  have  also  been  investigated  thoroughly  by  Lee 
(1993).  The  three  thermal  effects:  (i)  global  expansion  by  uniform  plate  temperature,  (ii) 
local  expansion  by  temperature  variation  over  the  plate,  and  (iii)  thermal  moment  induced 
by  temperature  gradient  across  the  plate  thickness  are  included  in  the  investigation.  The 
single  mode  Galerkin  mediod  and  the  EL  (Roberts  and  Spanos,  1990)  technique  are 
used.  The  anal5tical  continuum  approaches  have  been  so  far  limited  to  uniform  or  linear 
temperature  distributions  and  to  beams  and  rectangular  plates  of  either  simply  supported 
or  clamped  edges. 

For  over  three  decades,  the  finite  element  method  has  been  the  predominant  method 
for  structural  mechanics.  However,  there  are  only  few  studies  on  nonlinear  random 
response  of  structures  using  the  finite  element  method.  Hwang  and  Pi  (1972)  have 
investigated  a  simply  supported  rectangular  isotropic  plate  subjected  to  rain-drop  type 
uniform  intensity  random  acoustic  loads.  The  high  precision  18  degree-of-ffeedom  (DOF) 
triangular  plate  bending  element  developed  by  Cowper  et  al.  was  used.  Both  first  and 
second-order  nonlinear  stiffness  matrices  are  developed  to  account  for  large  deflections. 
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The  finite  element  nonlinear  equations  of  motion  were  treated  as  a  linearized  eigenvalue 
problem  with  an  iterative  scheme.  The  acceleration  spectra  at  the  plate  center  were 
obtained  at  three  pressure  levels.  No  comparison  was  made  with  other  approximate 
solutions. 

Busby  and  Weingarton  (1973)  used  the  finite  element  method  only  to  obtain  the 
nonlinear  differential  equations  of  motion  which  are  expressed  in  terms  of  the  normal 
mode  coordinates.  The  EL  method  is  then  used  for  the  solution  of  these  equations.  Mean 
square  deflections  at  the  midspan  are  obtained  for  beams  with  both  ends  simply  supported 
and  both  ends  clamped.  However,  comparisons  with  other  solutions  were  not  made. 

Chiang  (1988)  has  presented  a  finite  element  method  for  large  deflection  random 
response  of  beams,  plates  and  built-up  panels  subjected  to  acoustic  loads.  Geometrical 
stiffness  matrices  to  account  for  the  induced  inplane  forces  due  to  large  deflections  were 
developed  for  an  isotropic  beam  element  and  an  isotropic  rectangular  plate  element.  Root- 
mean-square  (RMS)  maximum  deflections  and  RMS  maximum  strains  are  obtained  for 
beams  and  rectangular  plates  with  simply  supported  and  clamped  boundary  conditions. 
The  finite  element  results  are  in  good  agreement  with  single-mode  Fokker-Planck- 
Kolmogorov  (FPK)  equation  and  analytical  equivalent  linearization  solutions. 

Locke  (1988)  and  Locke  and  Mei  (1990)  extended  the  finite  element  method  to 
isotropic  beam  and  rectangular  plate  structures  subjected  to  thermal  and  acoustic  loads 
applied  in  sequence.  The  thermal  load  considered  is  a  steady-state  temperature  distribu¬ 
tion  AT(x,y).  The  thermal  postbuckling  structural  problem  is  solved  first  to  obtain  the 
deflection  and  thermal  stresses.  The  deflection  and  thermal  stresses  are  then  treated  as 
initial  deflection  and  initial  stresses  in  the  subsequent  random  vibration  analysis.  The 
Newton-Raphson  iterative  mefliod  is  used  in  the  thermal  postbuckling  analysis.  For  the 
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nonlinear  random  vibration,  the  linear  mode  shapes  of  the  thermally  buckled  structure  are 
used  to  reduce  the  order  of  the  system  equations  of  motion  to  a  set  of  nonlinear  modal 
equations  of  a  much  smaller  order.  The  EL  technique  is  then  used  to  iteratively  obtain 
for  RMS  responses.  Excellent  agreement  has  been  obtained  between  the  finite  element 
and  the  Galerkin/numerical  simulation  results  by  Seide  and  Adami. 

Jay  Robison  (1990,  1991)  has  derived  a  numerical  integration  routine  from  a  set  of 
unified  single  step  integration  algorithms  using  a  weighted  satisfaction  of  the  equilibrium 
equation  governing  the  large  deflection  random  response  of  laminated  composite  plates. 
The  equilibrium  equations  are  derived  using  large  deflection  finite  element  formulations. 
The  in-plane  inertia  terms  are  considered  in  the  formulation,  however,  rotary  inertia 
terms  are  assumed  negligible.  Probability  density,  spectral  density  and  autocorrelation 
functions  of  the  maximum  displacement  and  strain  responses  are  presented  for  three 
acoustic  excitation  levels.  Classical  thin  plate  boundary  conditions  and  pseudo  white 
noise  excitation  are  used  in  this  investigation. 

Chen  (1990)  and  Chen  and  Yang  (1991)  have  presented  a  finite  element  formulation 
combined  with  stochastic  linearization  and  normal  mode  methods,  including  the  geo¬ 
metrical  nonlinearity  for  the  study  of  random  vibration  responses  of  beams,  frames  and 
composite  plates  subjected  to  simultaneously  spatial  and  temporal  Gaussian  stationary 
nonwhite  and  nonzero  mean  random  excitations. 

Chen  and  Mei(1993)  presented  a  finite  element  formulation,  solution  procedure  and 
results  of  a  study  attempted  to  analyze  nonlinear  random  response  of  beams  subjected  to 
acoustic  and  thermal  loads  applied  simultaneously. 

13  Outline  of  the  study 

The  acoustic  fatigue  life  prediction  problem  is  very  important  to  military  and  civil 
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aircraft.  The  problem  consists  of  three  major  parts:  (1)  acoustic  loading  analysis;  (2) 
determination  of  the  response  of  structures;  and  (3)  estimation  of  the  fatigue  life  of  the 
materials.  In  this  study,  only  the  second  part  is  concerned. 

Since  the  structural  response  due  to  intense  acoustic  pressure  levels  is  extremely 
nonlinear  and  strongly  dependent  on  the  thermal  environment,  the  method  of  superposition 
is  not  applicable.  Rather  the  thermal  effects  must  be  integrated  and  coupled  directly  into 
the  acoustic-structural  analysis.  The  finite  element  formulations  by  Chiang,  C.  T.  Chen 
and  Robinson  did  not  consider  temperature  effects.  The  formulation  by  Locke  treated  the 
two  loads  in  sequence;  thus  there  is  no  inter-dependence  between  the  thermal  effects  and 
the  acoustic-structural  response.  In  the  aforementioned  literature  survey,  it  appears  that 
studies  of  the  nonlinear  random  response  of  structures  subjected  to  simultaneously  applied 
acoustic  and  thermal  loads  using  finite  element  method  are  not  available  in  the  literature. 

In  addition,  only  normal  incidence  acoustic  pressure  loads  have  been  considered  in 
all  the  existing  investigations.  It  appears  that  studies  of  the  nonlinear  random  response 
of  plates  subjected  to  a  grazing  incidence  acoustic  wave  using  finite  element  method  are 
not  available  in  the  literature. 

Therefore,  this  dissertation  will  develop  a  finite  element  formulation  and  solution 
procedure  which  is  believed  to  be  the  first  attempt  to  analyze  nonlinear  random  response 
of  complex  composite  panels  subjected  to  simultaneous  acoustic  and  thermal  loads,  and 
the  acoustic  pressure  can  be  either  normal  incidence  or  grazing  incidence. 

In  Chapter  2,  the  formulation  for  the  problem  is  derived.  The  formulation  is  based 
on  von  Karman  large  deflection  theory  and  the  first  order  shear  deformation  theory.  In 
Chapter  3,  the  solution  procedure  is  described.  Using  linear  vibration  modes  of  the 
thermally  buckled  structure,  the  governing  equations  are  reduced  to  a  set  of  nonlinear 
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coupled  modal  equations.  The  equivalent  linearization  technique  is  employed,  because 
satisfactory  results  have  been  obtained  using  this  technique.  Finally,  in  order  to  uncouple 
the  linearized  modal  equations  of  motion,  the  modal  transformation  is  utilized  once  more. 
In  Chapter  4,  the  numerical  results  are  presented.  It  includes  thermal  buckling  and 
postbuckling  results,  random  response  and  acoustic-thermal  combined  response.  Chapter 
5  is  concerned  with  conclusion  and  further  work  suggestions. 
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Chapter  2 


FORMULATIONS 


In  this  chapter,  tiie  governing  nonlinear  equations  of  motion  are  derived  for  a  plate 
of  arbitrary  shape  subjected  to  a  set  of  simultaneously  applied  thermal  and  acoustic 
loads.  The  thermal  load  is  taken  to  be  an  arbitrary  distribution,  but  steady-state,  i.e., 
AT  =  AT(a;,j/,2r).  The  acoustic  loading  is  considered  to  be  a  stationary  Gaussian 
pressure  wave,  with  the  extension  such  that  in  order  to  include  in  the  travelling  wave 
a  combination  of  wavelengths,  the  pressure  at  any  one  point  is  random  but  the  whole 
pattern  moves  in  the  direction  A  (  normal  incidence  A=0°  and  grazing  incidence  A=90°) 
with  wave  traveling  speed  c. 

The  following  features  are  considered  in  the  formulation: 

(a) .  Initial  imperfection  deflection  Wo(x,y), 

(b) .  In-plane  initial  forces  {No}, 

(c) .  Arbitrary  temperature  distribution  AT(x,y,z), 

(d) .  Large  deflections  in  von  Karman  sense, 

(e) .  Composite  materials  with  transverse  shear  deformation,  and 

(f) .  Acoustic  waves  directed  with  an  inclination  angle  A. 

The  assumption  regarding  temperature  independent  material  properties  is  utilized  in 
this  study. 
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The  three-node  triangular  Mindlin  (MINS)  plate  element  with  improved  transverse 
shear  is  extended  and  employed  in  this  study.  The  element  was  initially  developed  by 
Tessler  and  Hughes  (1985).  This  simple  plate  element  of  five  degrees  of  freedom  per 
node  in  large-scale  finite  element  structural  analysis  and,  especially,  in  nonlinear  analysis, 
has  great  computational  advantage.  Tessler  and  Hughes  have  found  that  the  transverse 
shear  energy  was  the  major  cause  of  difficulty,  therefore  a  special  interpolation  scheme, 
anisoparametric  interpolation,  was  devised.  In  addition,  the  element  transverse  shear 
energy  was  further  enhanced  by  a  suitable  (element  appropriate)  shear  correction  factor. 
Based  on  extensive  numerical  testing,  MINS  is  an  excellent  element.  They  concluded: 
“Due  to  its  reliability,  economy,  and  good  stress  recovery,  it  may  be  regarded  as  a 
viable  candidate  for  extension  to  shell,  laminated  composite  and  nonlinear  analyses.” 
The  finite  element  formulation  is  described  as  follows  for  the  nonlinear  random  response 
of  composite  panels. 


2.1  Element  Displacement  Functions 

A  typical  triangular  plate  element  is  shown  in  Figs.  2.1  and  2.2  to  discretize  a 
rectangular  panel.  The  element  displacement  functions  used  in  the  derivation  of  the 
equations  of  motion  are: 

Ux  =  u{x,y,t)  +  z‘il;y{x,y,t) 

Uy  =  v{x,y,t)  +  zipxixjyit)  (2.1) 

Uz  -  w{x,y,t) 

where  Ux,Uy,Uz  are  the  three  displacement  components  at  any  point  in  the  element;  u, 
V,  w  are  the  displacements  of  the  middle  surface;  and  ipx  and  ipy  are  the  rotations  of  the 
normal  around  the  x  and  y  axes  due  to  bending  only. 
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The  node  displacement  vector  is  defined  as  follows: 


=  [[wi,W2,  ,  LV’xl ,  V’x2 ,  V’x3 ,  tpyl ,  0j/2 ,  ,  (2.2) 

lui,U2,U3,Vi,V2,V3\\ 

The  interpolation  functions  for  the  MINS  element  are 
w{x,y,t)  = 

(2.3) 

=  L^l,i^2,i3,Mi,M2,Af3jW 

{V’}  =  L^i?^2)^33  OjO,  0j{^}  (2.4) 

rpyix,  y,  t)  =  [H^^\ {V’}  =  [0, 0, 0,  ^1,6, 6J  W  (2.5) 

u{x,y,t)  =  ^23^330, 0,0j{iyTO}  (2.6) 


v{x,y,t)  =  =  L0,0,0,63  63  6J{it'm}  (2.7) 


where  ^25  ^3  are  the  area  coordinates,  and  the  transformation  between  x,y  and  is 


['  ■  ■if 

=  Xi  X2  X3  <  6 

.yi  y2  yz  J 

Uz 

Y  p2?/3  -  xzy2  y2  -  yz  xz-  X2' 
—  xzy\  -  xiyz  yz  -yi  xi  -  xz 
lxiy2  -  X2yi  yi  -  j/2  0:2  —  a:i . 


(2.8) 
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where  2A  =  (x2  -  xi)(j/3  -  J/i)  -  (»3  -  xi)(y2  -  J/i),  A  is  the  area  of  the  triangular 
element,  (xi,yi)  is  the  coordinates  of  the  node  i;  and 

Li  =  ^(63^4  -  &2-A^6),  L2  =  ^{hiN5  -  63A4) 

O  O 

Lz  =  ^(&2iV6  -  ftiiVs),  Ml  =  i(a2iV6  -  asN^) 

M2  =  ^(03-^4  -  aiNs),  Ms  =  ^{aiNs  -  a2Ns) 

O  O 

■^4  =  4^1^25  -^5  =  4^2 ^3,  A’e  =  4^3^i 


01=2:32,  02  =  Xi3,  03  =  X21, 


&1  =  2/23)  h  =  2/31)  ^3  =  2/12 


Xjj  —  3^2  Xj^  yij  yi  yj 


jee 


i^dA  =  2A 


k\l\m\ 

(2  ■]-  “(■  /  “h  ?7i)! 


(2.10) 


2.2  Nonlinear  Strain-Displacement  Relations 


The  von  Karman  strain-displacement  relations  are  given  as 


iTxy  ) 

where  {e° }  is  the  in-plane  strain  vector,  and  {k  }  is  the  curvature  vector  such  that 


U,x  I  1  Wl-x  I  I  W,xWo,x 

I  "^2]  i^i 

^,y  +  v,x  J  I  2w,xW,y  J  [  W,xWo,y  +  w,y Wo,x 


W,yWo,y 


(2.12) 


=  {^m}  +  {4}  +  {^0} 

f  l/’y.x 

{«}  =  \  V’x.y  >  (2.13) 

I  V’y.y  +  2^x,x  , 

where  subscripts  m,  b  and  0  denote  that  the  inplane  strain  components  are  due  to 


membrane,  bending  and  initial  imperfection  deflection  respectively.  The  shear  strain- 
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displacement  relations  are  given  by 


where  the  subscript  denotes  the  derivative. 

23  Constitutive  Law 


(2.14) 


For  the  k-th  layer  of  an  orthotropic  material  with  an  orientation  (j),  the  stress-strain 


relations  are 


1  Qii  Qi2  Qi6 

{cr}k  =  <  CTy  >  =  Qi2  Q22  Q26  ({e}  -  {'^ATlk)  +  {<^No} 

.  J  It  lQi6  Q26  Qee  J  k 

=  [<5]ifc({e}  -  {eArjjt)  +  too} 


(2.15) 


J  k 


Qu  Q45 
Q45  Q55 


=  [0.]ifc{7} 


where  the  ffee-expansion  thermal  strain  vector  is 


'ai 

— 2c5 

2cs  <  Q2 

lo 

c  =  cos  <f>,  s  =  sin  (j) 

where  [Q]j^  is  the  transformed  reduced  stiffness  matrix  for  the  k-th  lamina,  and  {ctno} 
is  the  initial  stress  vector  corresponding  to  {Nq}. 
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2.4  Resultant  Laminate  Forces  and  Moments 


The  resultant  forces,  moments  and  shear  forces  per  unit  length  acting  on  a  laminate 
are  obtained  by  integration  of  the  stresses  in  each  layer  through  the  laminate  thickness 


({N},  {M})  =  {C7}],(l,z)d2 


(2.17) 


and 


lyz 


dz 


(2.18) 


-h.  \ 
2 


J  k 


where  {N}  is  the  resultant  force,  {M>  is  the  moment  and  <R}  is  the  shear  force  vector. 
The  laminae  shear  stiffness  is 


1^.1  =  /  [t 

2 


h 

\  t-  - 

<544  <545 


<554  <555 


dz 


(2.19) 


For  a  lamina  in  the  material  axes,  the  so-called  reduced  stiffnesses  are 

El 


Qii  = 


<5i2  = 


1  -  1/121^21 
U12E2  V21E1 


1  —  1^121^21  1  —  ^121^21 
Eo 

Q22  =  -, - ^ - 

1  -  1/12  ^'21 

<566  =  <^12 


(2.20) 


Qu  =  (?23 
<555  =  <J13 

where  C?23)<ji3  are  the  shear  modulus  of  the  materials.  Then  the  transformed  reduced 
stiffnesses  are 


r  c2 


[0]  =  m-'Kim- 

2sc 
—2sc 


i-T 


c2 

SC  SC  c^  — 


[r)  = 

s  =  sin  y>,  c  —  cos 


(2.21) 
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The  above  equations  can  be  written  as 


N1  [A  B]  re° 

Ml  B  D  1  K 


J  J  [  K  J  1^  Mat 
{R}  =  [A.]{7} 


m 


(2.22) 


where  [A],  [B]  and  [D]  are  the  laminate  extensional,  extension-bending  and  bending 
stifftiess,  respectively,  and  [As]  is  the  laminate  shear  stiffness.  And  the  free-expansion 
thermal  resultant  force  and  moment  vectors  are 


2 

{{Nat},  Wat})  =  j  [Q]WAT}W,z)dz 


(2.23) 


and  {No}  is  a  known  initial  force  vector. 


2.5  The  Principle  of  Virtual  Work 


The  virtual  works  done  by  internal  and  external  forces  are 


<5Wint  =  -f  {^«}'^{M}  -h  a{67}'^{R})  dA  (2.24) 


•^Wext  =  /  [(!)w(p(x,y,t;  A)  -  phw,tt) 
Ja 


(2.25) 


-f  du(-/9hu  tt)  +  ^v(-/)h  V  tt)]dA 

where  p(x,y,t;A)  is  the  acoustic  excitation,  A  is  an  incidence  angle  for  normal  or  grazing 
acoustic  wave,  p  is  the  mass  density  of  the  laminate,  and  o;  is  a  shear  correction  factor, 


(2.26) 


1  0.5- 


It  is  assumed  that  the  rotatory  inertia  effect  is  neglected  for  relative  thin  plate  (a/h>50). 


Some  geometric  matrices  are  given  below: 


[0]  = 


{G}  = 


'w,x 

0  ■ 

Wolx 

0 

0 

W,y 

o 

11 

0 

'^ojy 

1 - 

g 

U>,X_ 

y^Oiy 

'^0'}X  _ 

w 


[C.^h]{wh}  + 


w 


{4)  =  {«4}  =  5«(1«1{G})  =  [«1{«G} 


(2.27) 


\Cm\  — 


J_ 


y2z 

yzi 

yi2 

0 

0 

0 

0 

0 

0 

®32 

a:i3 

a:2i 

(2.28) 

XZ2 

a:i3 

X21 

1/23 

yzi 

yi2. 

[C'^i]  = 


1 

2A 


y2z  yzi  yi2 

XZ2  Xl3  X21 


dx 


It? 


(2.29) 


1 


^‘tl}ipl2  ^il)il)lZ  ^^^15  CipijjiQ 


2A  ^■^^22  ^tptl^ZZ  ^^^24  ^ij>ijy2Z  (^-0^26  J 

(2.30) 

where  the  12  C^^y’s  are  related  to  the  area  coordinates  and  the  coordinates  of  the  three 
nodes,  and  they  are  listed  in  the  Appendix  A. 


[C'5]  = 


1 


0 

0 

0 

2/23 

2/31 

2/12' 

a:32 

a:i3 

a:2l 

0 

0 

0 

(2.31) 

.2/23 

2/31 

2/12 

a:32 

a:i3 

0:21 . 

[^7^1  = 


1 

®32  a:i3  2:21 

1 _ 

“  2A 

.2/23  2/31  2/12  . 

(2.32) 
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[^7^]  “  o 

_  C'^^21  +  ^1  C,I)^22  +  ^2  C'V’^23  +  ^3  <^^^24  C'V’^25  C'^^26 

C'^^12  ^*^^13  C'^^14  +  ^1  +  ^2  C'V’^16  +  ^3 

(2.33) 

Substituting  Eqs.  (2.28)-(2.33)  into  Eqs.  (2.12)-(2,14),  one  obtains 

{^e»f  =  {Sw„flC^f  +  {Sv>if{C^tf[ef  + 

(2.34) 

+{6tos}^[Ci,j]^[e/  + 

{iV}  =  Mle”}  +  !BJ{4  -  {JVat}  +  {«■.) 

=  M[C„J{to„}  +  i[A]P]([C*j]  {«,»}  +  [C^]  {V-})  +  Ml«]  [C^i]  {w,J  (2.35) 

+M[«1  [C«]  m  +  [B][C»]W}  -  {iVir}  +  {(V.} 

{S^f  =  {6i>f\Cif  (2.36) 

{M}  =  +  [0]W  -  {Maj.} 

=  [B](C„]{u.„}  +  i(B)[«]([C*t]{«,»}  +  [C*«]{V>})  +  [Bm[Ci^]{wi,}  (2.37) 

+[B|[«1[C#*]{V-.}  +  [B]lCi){V-}  -  {JWat) 

{«7}^  =  {S«hf[Cylf  +  (2.38) 

[R]  =  [A,]  [C,J  {»(}  +  [^]  [C,*]  {f}  (2.39) 

The  virtual  work  principle  gives 

S^Vint  —  SWext  (2.40) 

Finally,  from  Eq.  (2.24): 
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SWin,  =  J  ({«£”f{Ar}  +  {5/cf{Jlf}+a{«7f{i?))<iA 

A 

=  J [C^b]{'i^b}  +  2 [^0^]  {^} 

A 

+[4[«1  [c#]  {«»»<.}  +  I-4II«1  [Ch]  {<P>}  +  miCi]W  -  {Nat]  +  {JV.}) 
+{^W6}^[C'^6]^[5]^^[A][Cm]'{win}  +  ^[■^][^] 

+ W[«l  [c*s]  {«>».}  +  Ml«l  [C#*]  {M  +  [B\[Ci\{i>}  -  {JVat}  +  {N,}) 

+{#}’’ [C*^]  ([A][CJ{®„}  +  [c#]  {«>»}  + 5  WM  [Cw]  {0} 

+Ml«l  M  {»».}  +  [A][e]  [C**]  {^,}  +  [BIICjKV}  -  {JVat}  +  {No}) 
+{6wti}^  [C^bl^  [^of^  ([^][^m]{y^m}  +  2["^][^1[^#]{^6}  +  [^H’]{^} 

[C^i]  {“>».}  +  MM  [c^,^,]  {M  +  |Bll(?il{0}  -  {Nat)  +  {K}) 

+MM  +  MM  [C'^v]  {^0)  [BjK’bKi^]  —  {Nat]  +  {Nc]) 

+{#}^[CJ^([Bl[C„){w„}  +  i[B]M[Cw]{«>s}  +  Jl-BIM  [<?**]  {!/>} 

+1B]M  [c#]  {"'So}  +  [-BIM  [CiSis]  {M  +  l^>llC6l{V’}  -  {Mat]) 

[C,s]^  ([A.]  [C,s]  {w,}  +  lAol  [C,*]  {V-}) 

+a{<0)’'[C,J^([Ao]  [C,^]  {Wi)  +  (Aol  [C,*]  {MdA 

(2.41) 
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and  from  Eq.  (2.25): 


^  (2.42) 

(^p{x,y,t-,X)  -  ph([H^\{wb}  +  L-ffw0j{^})) 

-{Swjn}^{{Hu}[Hu\{wm})  “  {Swm}^{{Hv}[H^\{wm})^dA 

2.6  Element  Equations  of  Motion 

The  application  of  the  principle  of  virtual  work  to  derive  the  element  equations  of 
motion  and  the  element  matrices  is  lengthy  and  tedious.  There  are  total  56  terms  in  the 
expressions  of  virtual  work.  The  equations  of  motion  for  the  MINS  plate  element  can 
be  written  in  the  matrix  form  as 
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/ 

'0  0 

0 

0 

+ 

[^o]^b 

V 

.0 

[*:L  J 

}.^o\mh 

(’  [^ATAtIj  [kNAT]hijj  0 
[%at]^6  [%Ar]^  0 

0  0  0 


'Wb  > 

[^olbm 

\ 

[h]^ 

0  > 

0 

/ 

k  Wrri  J 

’  [%o]j 

[^Nolb^ 

o' 

\ 

[kNo]^b 

[%o]^ 

0 

0 

0 

0 

/ 

'Wb  ' 

0  [nl]j^  [nl]j^  ■ 

■[nlolj  [nlo]bi,  O' 

+ 

[«lo]^6  [«lo]^  0 

if;  >. 

o 

o 

o 

/ 

.■UJm  > 

'Wb  '' 

l/ 

+  2  ( 

+ 

[nlNb]rl,b  hl^■6]^  0 

0  0  o_ 

0  0  0 

/ 

'Wb  ' 

'Wb  ' 

'[n2\  [n2]j^  O' 

[n%b  0 

< 

0 

►  +  a 

[^s]i  [^s]i^  0 

[^s]^b  0 

-! 

0  ^ 

.  0  0  0_ 

>.  '^m  > 

.  0  0  0_ 

'Wb  ' 

Hfc 

0 

+ 

0 

0 

0 

0 

Nm. 

fO  > 

'{PATo}b  ' 

fO  >1 

'{PNoo}b  " 

►  +  < 

{pat}^ 

>  +  ^ 

{PATo}.^ 

^  ^ 

0 

►  +  •< 

{PNoo}^ 

.0 

{PAT  }m 

.0 

^  {PNo}m  ‘‘ 

.0 

(2.43) 
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or  in  the  short  form, 

([k]  +  [ko]  -  [kNAT]  +  [kNo]){w} 

+  +  [nlo]  +  [nlNm]  +  [nlNb]){w} 

I  (2.44) 

+-[n2]{w}  +  Q![ks]{w}  +  [m]{w} 

=  {Pp(*)}  +  {pat}  +  (paTo)  +  {PNo}  +  {PNoo} 
where  [m],  Dc]  and  {p}  denote  the  element  mass,  linear  stiffness  matrices  and  load 
vector,  respectively,  and  [nl]  and  [n2]  denote  the  first  and  second-order  nonlinear 
stiffness  matrices,  respectively.  The  subscripts  b,  ^  and  m  denote  the  transverse, 
rotation  and  in-plane  components,  respectively.  The  subscripts  s,  o.  No,  NAT,  Nm, 
Nb  denote  the  stiffness  matrices  which  are  due  to  transverse  shear,  Wo(x,y),  {No}, 
■(Nin}(=[A]{e^»,  and  {Nb)(=[B]{/c}),  respectively.  The  expressions  for  element 
stiffness,  nonlinear  stiffness,  mass  matrices  and  load  vectors  are  given  in  the  Appendix  B. 

2.7  Grazing  Acoustic  Wave 

The  grazing  incidence  acoustic  wave  applied  on  a  plate  can  be  treated  as  an  extension 
of  the  plane  wave  which  includes  in  the  travelling  wave  a  combination  of  wavelengths 
such  that  the  pressure  at  any  one  point  is  random  but  the  whole  pattern  moves  in  the 
direction  A  with  wave  traveling  speed  c.  The  pressure  distribution  on  a  plate  is  then 
given  by  [Clarkson] 


CC 

p{x,y,t)  =  j-  j  (2.45) 

— oo 

where  x  is  the  coordinate  along  the  wave  travelling  direction,  and  assume  that  the 
pressure  distribution  is  independent  of  y.  This  model  is  suitable  to  represent  the  waves 
in  the  progressive-wave  test  facility. 
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For  a  random  plane  wave,  the  pressure  at  a  point  can  be  written  in  the  Fourier 
integral  form; 

00 

p(t)  =  ^  y  P{u)e^‘^^dw  (2.46) 

— OO 

Physically,  this  is  equivalent  to  representing  the  oncoming  wave  as  the  sum  of  an  infinite 
number  of  waves  of  different  wavelength. 


The  nodal  force  vector  of  an  element  can  be  calculated  as: 

{Pp{t)}  =  J  {Rw}p(x,y,t‘,\)dA 
A 

00 

=  ±  y  P(u)e‘“‘  J  c-^^>'{n„}dAdu 


{!"(")}  =  / 


e  -  {B.w}dA 


then  we  have 


OO 

{Pp(i)}  =  ^  /  P(‘o){Y(ic)]e^*du, 


(2.47) 


(2.48) 


(2.49) 


where  {Ht„}  is  defined  in  Eq.  (2.3)  as 

^  r 

w(x,y,t)  =  ^  =  {H^ix,y)f{wi}  +  {H^^]  W  (2.50) 

W  J 

and  the  nodal  displacement  components  {wj,}  and  {tp}  are  defined  in  Eq.  (2.2).  Therefore, 
P(a;){y(a;)}  is  the  Fourier  transform  of  {pp(<)}.  The  spectmm  density  of  {pp(f)}  is 
(Clarkson,  1986) 

{SX")}  =  lim 

2^^°°  T  (2.51) 


25 


Fig.  2.4  Grazing  Acoustic  Wave 


Chapter  3 


SOLUTION  PROCEDURE 

3.1  Static  Component  and  Dynamic  Component 

The  system  equation  of  motion  can  be  written  as 

[M]  {ty}  +  ([if]  +  [K.]  -  [^fjVATl  +  [ifAro]){»'} 

+i[JVl]{»'}  +  i[iV2){»y}  +  a[K,]{W}  (3.1) 

=  {Ppit)}  +  {-Pat}  +  {Pato}  +  {Pno}  +  {Pnoo} 
where  { =  [Wi,  Wm\  •  The  subscripts  s,  o.  No,  NAT,  Nm,  Nb  denote  the  stiffness 

matrices  which  are  due  to  transverse  shear,  Wo(x,y),  {No},  {N^^>,  {Nin>(=[A]{e^», 

and  {N5>(=[B]{/c»,  respectively.  And 

[^1]=  (M]  +  +  Mjvj]  +  MoD  (3.2) 

assembly 

In  order  to  solve  the  system  equation  of  motion  Eq.  (3.1),  an  innovative  solution 

procedure  is  described  as  follows.  First,  the  response  is  assumed  to  be  the  sum  of 

and  {W}^,  i.e.  {W}  =  {W}^  +  {W}^,  where  {W}^  denotes  the  time-independent  or 

static  component  and  { W}^  the  time-dependent  or  dynamic  component.  The  displacement 

represents  a  stable  static  equilibrium  position  due  to  thermal  load  and  the  mean 

of  the  random  excitation  E[{Pp(i)}];  while  {1^}^  represents  the  zero-mean  random 

response.  Substituting  {W}  into  Eq.  (3.1)  and  regrouping  the  terms,  the  equation  is 

of  the  form  ^’({H^}^)  =  G{{W}^)  with 

-P({^}s)  =  ([-^]  +  o\Ks]  +  [Ko]  -  [iif^■A^]  +  [A'jvo]){^}s 

+5[-'Vl],{»'}.  +  i[Ar21„{lV}.  -  {Pat}  -  £KP,(i)}) 
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aim,) = iM]{w]^ + ((/<■] + a[K,\ + m  -  i/rjfiT) + [KNc])m, 

+|[M),{W^},  +  j(lJV21„  +  (Af21,.  +  [JV2]„)({M'},  +  {W},) 

-{PS)) + E[{p^m 

In  Eq.  (3.3),  the  load  vectors  {P^To},  {Pno}  and  {PjVoo}  are  temporarily  dropped.  The 
left-hand-side  of  the  equation,  is  independent  of  time  t;  while  the  right-hand- 

side,  is  time  dependent.  Therefore,  the  only  possibility  for  both  F  and  G  to 

exist  is  that  both  F  and  G  equal  to  zero,  and  the  following  two  equations  are  thus  obtained: 
m  -f  a[K,]  -F  [Ka]  -  [Knat]  + 

+5[M].{W'}.  +  j!^2]„{W'}.  =  {i’AT}  +  £:l{P,(i)}) 

lM][w}^  +  ([A-j  +  q[a:,)  +  [ir,l  -  [A-jvat]  +  \Ki)S{V), 

+  ( {V),  +  5[iVl],{H'}, 

+51^11, {W'},  +  5(IJV2]„  +  [JV2]„  +  [JV21„)({ir},  +  {W},) 

=  {PS)}  -  E[{PS)]\ 

where  the  subscript  [  denotes  tiiat  the  corresponding  stiffness  matrix  is  evaluated  with 
the  static  deflection,  and  [  is  evaluated  with  the  dynamic  deflection. 

Examination  of  Eqs.  (3.5)  and  (3.6)  reveals  that  both  equations  are  nonlinear.  The 
Newton-Raphson  iterative  method  is  used  to  determine  {W}^  from  Eq.  (3.5).  Combined 

T 

normal  mode  method  and  equivalent  linearization  technique  are  then  applied  to  Eq.  (3.6) 
to  obtain  {W}^. 

In  Eq.  (3.6),  the  matrix  [Nl]s  which  is  ignored  in  the  reference  (Locke,  1988)  is  due 


to  the  membrane  component  of  thermally  postbuckling  displacement  {Wat}.  The  term 
[N2]ss{W}t  is  due  to  the  combined  effect  of  {Wat}  and  {Wat},  which  has  a  different 
coefficient  in  Locke’s  formulation(1988). 


When  comparing  with  the  Sequential  LazA  method(SQL)(Locke,  1988),  this 
SiMultaneous  Load  method(SML)  is  much  logical  mathematically  and  straightforward, 
it  is  then  easier  to  formulate  nonlinear  problems  with  combined  loading.  The  solution 
procedure  itself  can  take  care  of  the  inter-dependence  between  the  thermal  effects  and 
the  acoustic-structural  response.  In  the  SQL  method  engineering  judgment  is  essential, 
otherwise  some  terms  might  be  missed. 

3.2  Thermal  Buckling  and  Large  Thermal  Deflection 

To  obtain  the  critical  buckling  temperature  change  AT^r  (x,  y)  for  a  plate,  it  is  assumed 
that  the  prebuckling  configuration  is  flat  and  without  coupling  between  bending  and 
extension,  then  the  linear  system  equation  in  membrane 


[Km]{W^)  =  {P^Ar} 


(3.7) 


is  solved  for  an  assumed  temperature  distribution  AT  first.  Once  {Wm}  is  obtained,  one 
can  calculate  the  first-order  nonlinear  stiffness  matrix  [A^l^rm].  Because  the  nonlinear 
term  due  to  transverse  deflection  does  not  exist,  the  stability  equation  of  the  investigated 
plate  becomes 


where 


{[Ki]  -  [Knat]  +  [.ATI^Vto])  <  >  =  0 

[a^  J 


(3.8) 


0  0  ■ 

+  a 

.0 

(3.9) 


Examination  of  Eq.  (3.8)  reveals  that  [Ki]  is  independent  of  temperature  change,  and 
(-&At]  +  [NlNm])  is  proportional  to  temperature  change.  Therefore  Eq.  (3.8) 
describes  an  eigenvalue  problem.  The  critical  buckling  temperature  change  is 


ATcr  =  AiAT 


(3.10) 
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where  Ai  is  the  lowest  eigenvalue.  For  a  non-symmetrically  laminated  plate,  the  ATcr 
from  Eq.  (3.10)  is  referred  as  reference  temperature,  ATre/. 

The  iterative  solution  scheme  is  to  seek  a  solution  of  =  0.  The 

Newton-Raphson  iterative  method  is  a  well  established  procedure  for  solving  time- 
independent  nonlinear  problems.  This  method  involves  a  repeated  solution  of  the  equation 
for  the  i-th  iteration 


(3.11) 

Then  and  {AP},  are  updated  by  using  =  {W},^i  +  {AlV},_j^,j. 

The  solution  process  seeks  to  reduce  the  load  imbalance,  and  consequently  {AW}^,  to 


a  specified  small  quantity.  The  tangent-stiffness  matrix  is  determined  from 

■dF({vr},)' 


[K]t  — 


=  [K\  +  a[Ks\  +  [Ko 


(3.12) 


+  [Kno]  -  [Knat]  +  +  [-^2]^^ 

and  the  load  imbalance  vector  is 

{AP},  =  {P},  -  ([dft:]  +  a[Ks]  +  [Ko]  -I-  [Kno] 

-  IKnat]  +  jIlVl],  +  |[JV2]„){ir}. 
where  {P}^  =  {Pat}  +  ^[{-PpCO}]-  Th®  linear  buckling  mode  shape  from  Eq.(3.8) 

multiplying  a  scale  factor  is  usually  taken  to  be  the  initial  trial  solution  of  Eq.  (3.11). 


(3.13) 


33  Nonlinear  Random  Response 


After  solving  {W}^  from  Eq.  (3.5)  and  evaluating  the  matrices  [Nl]^  and 
Eq.  (3.6)  is  ready  to  be  solved.  Firstly,  Eq.  (3.6)  will  be  reduced  to  a  system  of  coupled 
nonlinear  modal  equations  with  reduced  degree-of-freedom.  The  linear  vibration  modes 
of  the  deformed  structure  are  used  to  transform  the  system  equation  of  motion  to  modal 
coordinates.  The  resulting  nonlinear  modal  equations  of  motion  are  then  linearized  using 
the  equivalent  linearization  method.  Finally,  in  order  to  uncouple  the  linearized  equations 
of  motion,  a  modal  transformation  is  used  once  more. 
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Coupled  nonlinear  modal  equations 


From  Eq.(3.6),  the  linear  frequencies  and  mode  shapes  of  the  deformed  structure  can 
be  obtained  by  solving  the  eigenvalue  problem 

([A-J  +  a[Jir.j  +  [K,\  -  [Kkat]  +  [iTjf,]  +  [ffl],  +  [JV2]„){^}„ 

where  it  is  assumed  that 


(3.14) 


+  IfiVllJW'},  (3.15) 

and 

+  j[JV2]„{»'},  +  ilJV2]„{H^},  (3.16) 

Actually,  from  numerical  test  this  assumption  doesn’t  introduce  significant  error.  Solving 
Eq.  (3.14),  the  tmncated  modal  matrix  is  given  by 


[^]  —  [{^}l5  {^}2j  •  •  •  5  (3.17) 

where  N  is  the  number  of  modes  to  be  used  for  the  analysis  of  Eq.  (3.6).  Now  {W}^ 
can  be  written  in  terms  of  the  modal  amplitudes  as 

N 

=  [<t>]{q]  =  {<i>}^qn  (3.18) 

71=1 

Using  Eq.  (3.18),  the  first-order  nonlinear  system  stiffness  matrix  can  be  written  in  terms 
of  {g}  as  the  sum  of  first-order  nonlinear  system  modal  stiffness  matrices 

[^l],  =  i;?n[iVllS">  (3.19) 

71  =  1 

and 

[Ari]<"'=  ^  ([nlj  +  [nl«„l  +  [nl,]  +  [nlOTl)<">  (3.20) 

asse.mbly 
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(3.21) 


j([A'2]„  +  [iV21„){W'},  =  \[K2\„{q)  (3.22) 

1  1  ^ 

■^muw},  =  -Y^  !„[/aiir*{s}  (3.23) 

n=l 

where  the  ith  column  in  [Kl]^,[K2]^^  and  [K2]^^  are 

mhi  =  [(Vljf>{H'}.;  {K2},„  =  ([)V2]„  +  [iV2],j(«{)r}, 

(3.24) 

{iir2}W  =  iiV2]<;>{H^}, 

The  second-order  nonlinear  system  stiffness  matrices  as  the  sum  of  second-order 
nonlinear  system  modal  stiffness  matrices  are 


!(V2].,  =  f;,„lAr2lW 


n=l 


and 


similarly 


[JV2jW=  Y. 

assembly 


=  E?4(V2]£> 

n=l 

[(V2]W=  Y 

assembly 


m\„  =  YYmr[N2t’'> 

n=l  r=l 

[(V2]<r'>=  E  [n2]<r> 

assembly 


(3.25) 

(3.26) 

(3.27) 

(3.28) 

(3.29) 

(3.30) 


where  superscripts  [  and  [  denote  the  corresponding  nonlinear  modal  stiffness 
matrix  is  evaluated  with  the  modes  and  {^}^. 
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Substituting  these  nonlinear  modal  stiffness  matrices  into  Eq.  (3.6),  a  set  of  nonlinear 
coupled  modal  equations  of  motion  can  be  expressed  as 

{4>fG({W},)  =  IMM  +  lKiun,.A<!} 

n=l 

+  +  i[JV2lW  +  i(JV2i(:))M{,} 

+  W{?}  -{/}=(! 

\  n=l r=l  / 

where  the  modal  force  vector  and  the  diagonal  modal  mass  and  linear  stifftiess  matrices 
are 


{/}  =  W{P}t 

mi 


fMJ  =  [<f>Y  [MM  = 


0 


mjv 


(3.32) 

(3.33) 


\^llinear  -  [(f>f  (^[K]  +  +  [Ko]  -  +  [Kno]  +  [<;6] 


L  0 


{p}t={Ppm-mPp{m 


(3.34) 

(3.35) 


and  [  J  denotes  a  diagonal  matrix.  Equation  (3.31),  is  a  set  of  coupled  nonlinear  modal 
equations,  it  is  then  linearized  by  using  the  equivalent  linearization  method  (Atalik  and 
Utku,  1976,  and  Roberts  and  Spanos,  1990). 

Equivalent  linearization  and  coupled  linearized  modal  equations 


Rewriting  Eq.  (3.31)  in  the  form 


{^({?})}  +  =  {/}  (336) 
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The  corresponding  linear  form  of  Eq.  (3.31)  can  be  expressed  as 


[£■]{,} +rMJ  {«}  =  {/}  (3.37) 

where  [.^]  is  an  equivalent  linear  modal  stiffness  matrix.  The  error  vector  involved 
in  using  Eq.  (3.37)  instead  of  Eq.  (3.36)  is  given  by  the  difference  between  the  two 
equations  as 

{e}  =  MW)}  -  [^H?}  (3.38) 

The  equivalent  linear  stiffness  matrix  [.^]  can  be  found  by  requiring  that  the  mean  square 
value  of  error  be  a  minimum,  thus  we  have 

ai;f{e}^{e}] 

- n,r  =  l,2,...,iV  (3.39) 

Substituting  {e}  and  {£(({9})}  into  Eq.  (3.39),  the  equivalent  stifihiess  matrix  [K\  can 
be  determined  from  the  equation 

-E[{9}{?}’'][^]  (3.40) 

The  right-hand  side  of  Eq.  (3.40)  can  be  evaluated  as 
E  [{9}{9}"']  =  E  [{,}{5}^]  \K\  +  E 

n 

r  /  T 

+  E'E[9"{9}<9}1  [l^l"'(j[(Vl]i“>  +  i[iV2]W  +  i[Ar21<:>j[al] 

n  r 

=  £[{4}{«}"']  +£[{9}{9}"']  +  5(*^2]„i]"' 

+  5EE^[m'W(9}’"] 

n  T 

(3.41) 
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since  E 


(ln{q}{q}^ 


=  0  for  Gaussian  process,  and 


(3.42) 

If  the  covariance  matrix  £  |{g}{9}^]  is  known,  the  equivalent  linear  stiffness  matrix  [R] 
can  be  determined  from  Eq.  (3.41).  However,  {g}  has  to  be  obtained  first  from  Eq.  (3.37) 
and  \K'\  is  not  known.  In  order  to  solve  Eq.  (3.37),  modal  coordinate  transform  is  used 
once  again  with  an  iterative  scheme. 

Uncoupled  linear  modal  equations 


The  modal  transformation  is  used  once  more  and  it  is  determined  from  the  equation 


(3.43) 

and  the  modal  transformation  matrix  is  defined  as 

(3.44) 

Equation  (3.37)  thus  becomes  a  set  of  uncoupled  modal  equations 

rij  +  CVj  +  njrjj  =  fj  j  =  l,2,...,N 

(3.45) 

where 

=  ^ — 

nrij 

(3.46) 

fjp 

(3.47) 

2nd  (lij  =  2(luiir]j  is  the  modal  damping  term  which  has  been  added  to  Eq.  (3.45),  and 
02/1  is  the  first  linear  frequency  of  the  deformed  system,  which  is  obtained  from  Eq.  (3.14). 


Solutions  for  the  uncoupled  modal  equation  of  motion  Eq.  (3.45)  for  the  case  of 
Gaussian  white  noise  uniform  random  load  p{t)  are  given  in  the  following  form 

E[vm]  =  Spfjfkijk  (3.48) 
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where 


00 


Att^ 


—00 


(3.49) 


(3.50) 


where 

Hjiuj)  =  —5 - ^ — — - 

fiy  —  +  l^UJ 

and  Sp  is  the  double  sided  loading  spectrum  density.  Using  Eq.  (3.44),  the  covariance 
matrix  of  the  coupled  modal  amplitude  becomes 

•-iTl 


(3.51) 


The  deflection  spectrum  density  is  as  follows: 


Gw  =  SpfjhHj(^)Ht(-w) 

-iT 


(3,52) 


For  grazing  incidence,  flie  covariance  terms  for  the  case  of  ideal  white  noise  excitation 
can  be  expressed  as 

CX) 

E[rjjrik]='2.  J  Sp{u)\Zj{u)Zk{c^)\\Hj{u;)Hki-uj)\duj  (3.53) 


where 


Zj{u)  = 

rrij  ^  J  j 


(3.54) 


and 


{yM}=  y:  /fco.(^^)-.sin( 

element''.  L  V  ^  /  \ 


tjxsin^' 


'Ey, 


■dA  (3.55) 


In  order  to  calculate  E[rijr}k]  from  Eq.  (3.53),  a  simple  numerical  integration  method 
is  used,  the  cut  off  frequency  is  1.5  *  fijv- 
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Iterative  solution  procedure 


Therefore,  Eqs.  (3.40),  (3.41)  and  (3.51)  can  be  used  to  determine  and 
However,  since  each  of  these  quantities  is  dependent  on  the  other  and 
these  equations  are  nonlinear,  consequently,  they  must  be  solved  by  an  iterative  method. 

The  first  approximation  of  [/?]  and  E  is  obtained  by  neglecting  the  cross 

terms  in  ,  (i.e.  E[qnqr]  =  0  for  n  ^  r),  and  assuming  all  the  equations  in 

Eq.(3.37)  are  uncoupled.  The  diagonal  terms  in  the  equivalent  linear  stiffness  matrix  [R] 
can  be  expressed  from  Eqs.  (3.40)  and  (3.41)  as 


Knn=Kr,  +  HK2nnf^'>E[ql] 


(3.56) 


where  {K2nn)^^^^  are  the  diagonal  terms  of  the  second-order  nonlinear  stiffness  matrix, 
the  subscript  nn  denotes  the  diagonal  term  and  the  superscript  (nn)  denotes  the  term  is 
due  to  the  n-th  mode  From  Eq.  (3.43)  can  be  written  as 

\2  ^nn 


nt 


TUr. 


(3.57) 


and  E  [§„]  can  be  found  from  Eq.  (3.48)  for  this  uncoupled  first  approximation  (77^  = 

fn  =  fn/rUn)  tO  be 


Using  Eqs.  (3.57)  and  (3.58),  ^£^[9^]  can  be  determined  to  be 

E[ql]  =  (4^2  +  4(0,  _  5^/2 


(3.58) 


(3.59) 


where 


C  =  Sp/„V/3^m.(/f2.,)(”") 


(3.60) 
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To  begin  the  iteration  process  for  the  cross  terms  are  no  longer  neglected,  the  cross 
correlation  terms  Elqnqr]  can  be  evaluated  using  Eq.  (3.48)  with  Eq.  (3.49)  and  (3.51). 
Using  Eqs.  (3.40)  and  (3.41),  the  equivalent  linear  stiffiiess  can  be  computed.  After 
obtaining  [^] ,  a  new  iterative  cycle  begins.  The  iterative  process  goes  on,  until  some 
convergence  criteria  are  satisfied,  —  ”  — ~  —  iteration  is 

terminated.  The  covariance  matrix  can  be  found  from 


(3.61) 
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3.4  Strain  Formulas 


After  the  displacements  for  a  given  combination  of  thermal  and  acoustic  load 
condition  are  known,  the  element  strains  can  be  calculated  using  Eqs.  (2.11),  (2.12), 
(2.13)  and  (2.14).  Because  the  displacements  consist  of  two  parts,  dynamic  and  static. 
The  strains  can  be  expressed  as 


where 


and 


where 


{€}  =  {e°}  +  2{«}  =  {e}3  +  {e}i 
{7}  =  {7)3  +  {l}t 

+  2  (  [^#]  {^}s) 

=  X!  +  X  X 

j=l  i=l  jt=l 


(3.62) 

(3.63) 


(3.64) 


(3.65) 


{42i*;  =  'o^^]i{[^i)h]{<i>h}k  +  -Ip-Ip] 


[^].  = 


4>hi,x  0 

0  4>bj,y 

<f>bj  iy  4*bj  tx 


For  the  shear  strain. 


(3.66) 

(3.67) 

(3.68) 


{7}s  =  +  [^‘iip]{i’]i 


(3.69) 


and 


N 

{7}*  =  X  ([*^76]  {^b]j  +  [Cyi]  =  Y,  {7}iy9i  (3.70) 

3=1 


N 

c 

j=l 
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E[qiqiqk<n]  =  E[qiqj]E[qkqi]  +  E[qiqk]E[qjqi]  +  E[qiqi]E[qjqk]  (3.73) 
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Chapter  4 


NUMERICAL  RESULTS  AND  DISCUSSION 

4.1  Formulation  and  Computer  Program  Validation 

4.1.1.  Twisted  Ribbon  Tests 

A  cantilevered,  thin  rectangular  plate  subject  to  a  twisting  moment  at  the  free 
end  is  regarded  as  a  severe  test  for  plate  bending  elements  under  large  aspect  ratio 
distortions(Robinson,  1979).  Herein,  the  mesh  A  (Figs.  4.1  and  4.2)  results  has  been 
repeated.  The  tip  deflection  results  are  compared  with  the  results  produced  by  the  popular 
nine  degrees-of-freedom,  thin  triangles,  namely,  BCIZl,  HSM,  HCT  and  DKT,  are  shown 
in  Figs.  4.1  and  4.2.  In  this  comparison,  the  MINS  which  is  adopted  for  the  present 
random  analysis  appears  to  outperform  the  other  triangular  elements. 

4.1.2.  Static  Analysis  of  a  Rhombic  Cantilever 

This  problem  deals  with  the  analysis  of  a  rhombic  cantilevered  plate  subjected  to  a 
uniform  load.  The  geometry  and  material  properties  are  given  in  Figure  4.3.  Experimental 
results  of  this  problem  are  available  for  comparison  (Clough  and  Tocher,  1965).  A  4  by 
4  mesh  (32  elements)  is  used  and  the  results  obtained  with  MINS  and  test  are  given  in 
Table  4.1.  It  is  observed  that  even  with  this  coarse  n»sh,  the  MINS  element  gives  results 
in  good  agreement  with  the  experimental  values  (Batoz  et  al.,  1980). 


42 


4.13.  Thermal  Buckling  and  Post  Bucklini 

For  this  problem,  a  square  plate  with  two  boundary  conditions  is  calculated,  one  is 
simply  supported  and  the  other  is  clamped.  The  results  are  compared  with  the  results  by 
Paul  (1982)  and  Singh  (1993),  the  agreement  is  very  good. 

For  the  simply  supported  square  plate,  the  results  are  given  below: 


(4.1) 


(4.2) 


The  plate  used  by  Chiang  (1988)  was  analyzed  as  a  validation  example.  The  plate 

is  of  the  following  dimension  and  material  properties; 

Young’s  modulus  E  =  10.5  x  10®psi 

mass  density  p  =  0.2588  x  W-Hhse(? jin^ 

damping  ratio  (  =  0.01 

length  G  =  15m 

width  b  =  12in 

thickness  h  =  0.04m 

Poisson’s  ratio  u  =  0.33 


The  results  obtained  by  the  present  study  and  Chiang  are  shown  in  Table  4.2.  Because 
the  elements  used  are  different,  the  difference  of  strain  values  for  N=4  is  relatively  large. 


The  others  are  close. 


4.2  Thermal  Buckling  and  Postbuckling  Results 

In  order  to  understand  the  response  of  a  plate  to  combined  acoustic  and  thermal  load 
better,  the  thermal  buckling  and  postbuckling  behavior  of  a  composite  plate  is  studied 
in  this  section.  Seven  cases  are  investigated.  The  graphite-epoxy  material  properties 
are  taken  as: 

Young’s  moduli  =  22.5  x  10®,  E2  =  1.17  x  10®ps^ 

Shear  moduli  G2Z  =  0.4  x  10®,  =  (?i3  =  0.66  x  10®p5i 

Poisson’s  ratio  =  0.22 

Therm,  expan.  coeff.  =  _o.04  x  10“®/f’°,Q:2  =  16.7  x  10~®/i^‘’ 

The  finite  element  results  are  presented  as  follows: 

4.2.1.  Effect  of  Extension  and  Bending  Couptinp 

The  bending  and  extension  coupling  is  studied  first.  A  two-layer  Gr/Ep  rectangular 
laminate  (  15x12x0.048  in.)  with  the  stacking  sequence  of  (0/90)  is  considered.  For 
this  case  the  extension  and  bending  coupling  matrix  [B]  is  not  equal  to  zero,  therefore 
the  critical  buckling  temperature  does  not  exist.  When  the  laminate  is  subjected  to 
the  temperature  change,  the  bending  deflection  occurred  immediately.  The  postbuckling 
deflection  is  shown  in  Fig.  4.4  for  the  simply  supported  boundary  condition.  The  AT^f 
used  in  the  figure  is  13.37F°  only  for  reference  purpose,  and  it  has  no  physical  meaning. 
In  the  calculation,  the  full  plate  is  discretised  to  128  elements  or  8x8x2  mesh.  From 
the  figure  one  can  see  that  there  is  no  bifurcation  critical  temperature. 
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42.2.  Isosceles  lyiangular  Plate 


The  second  problem  investigated  is  an  isosceles  right  triangular  plate  with  symmet¬ 
rical  stacking  sequence  of  (0/45/-45/90)s,  the  lengtii  of  two  perpendicular  sides  is  12 
in.  The  plate  is  simulated  by  144  elements.  In  Fig.  4.5,  the  postbuckling  behavior  is 
shown  for  simply  supported  and  clamped  boundary  conditions.  These  two  conditions  are 
theoretically  idealized,  the  boundary  conditions  in  the  real  world  are  somewhere  between 
them,  therefore  they  can  be  considered  as  the  upper  and  lower  bounds. 

42.3.  Effect  of  Shear  Deformation 

The  plates  studied  are  the  same  as  the  first  case,  i.e.  15x  12  in.,  but  the  stacking 
sequence  is  (0/45/-45/90)s.  The  ratio  of  length  to  thickness  investigated  are  a/h=312.5, 
200, 100,  50  and  20.  The  critical  buckling  temperature  is  shown  in  Fig.  4.6.  As  expected, 
the  results  show  that  when  a/h  is  greater  than  100,  the  shear  deformation  can  be  neglected. 
But  for  thick  laminates  the  shear  deformation  is  important. 

4.2.4,  Effect  of  Number  of  Layers 

The  dimension  of  the  plate  studied  in  this  problem  is  15x  12x0.08  in.,  and  the  mesh 
used  is  8x8x2  (128  elements).  The  boundary  support  condition  is  simply  supported. 
The  plate  consists  of  (45/-45)n.  In  Fig.  4.7,  it  is  shown  that  the  increase  of  number  of 
layers  reduces  the  response  due  to  the  reduction  of  the  extension  and  bending  coupling. 

4.2.5.  Postbuckling  Mode  Change 

The  fifth  problem  studied  is  a  36x12x0.048  in.  long  rectangular  laminate,  the 
stacking  sequence  is  (60/-60).  The  full  plate  is  modeled  by  a  18x6x2  mesh  i.e,  216 
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elements.  The  results  are  shown  in  Figs.  4.8  and  4.9.  The  laminate  is  subjected  only 
to  a  uniform  temperature  change  without  transverse  mechanical  load.  For  this  load  case, 
the  postbuckling  deflection  is  close  to  a  (3,1)  mode  shape  at  low  temperature.  When 
the  temperature  change  AT  is  greater  than  20°F,  there  is  a  mode  shape  change  and  the 
deflection  is  close  to  a  (4,1)  mode. 


4.2,6.  Thermo-Mechanical  Postbuckling 

The  plate  investigated  in  the  sixth  problem  has  the  same  dimension  and  stacking 
sequence  as  the  previous  problem.  The  mechanical  load  of  uniformly  distributed  0.01 
psi  is  applied  simultaneously  with  the  uniform  temperature  change.  Figures  4.10  and 
4.11  show  the  thermo-mechanical  deflection  and  the  deflection  shapes.  It  is  interesting 
to  note  that  the  maximum  deflection  exhibits  slightly  soften  behavior  at  low  temperature 
change  due  to  increasing  in  thermal  compressive  in-plane  forces.  The  deflection  shape 
at  this  low  temperature  is  a  combination  of  (1,1)  and  (3,1)  modes  due  to  the  presence  of 
mechanical  load.  However,  die  deflection  at  the  high  temperature  is  changed  to  a  (5,1) 
dominated  mode  shape.  The  mechanical  load  is  to  simulate  the  static  pressure  difference 
applied  to  the  aircraft  skin  panels. 

4A.7.  Effect  of  the  Skew  Angle  of  the  Plate 

In  this  problem,  the  skew  angle  /3  (see  Fig.  4.12)  of  the  plate  varies,  but  the  height 
(12  in.)  of  the  parallelogram  keeps  the  same.  The  length  of  the  plate  is  15  in.  and  the 
thickness  is  0.048  in.  The  height  of  the  plates  studied  is  equal  to  12  in.  Fig.  4.12  shows 
the  postbuckling  response.  When  the  skew  angle  /?  increases  the  deflection  reduces,  this 
is  due  to  that  the  length  of  90°  fibers  are  relatively  shorter  thus  making  the  plate  stiffen 
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4^.8.  Conclusions 


Seven  problems  were  studied  in  this  section.  The  extension  and  bending  coupling 
stiffness  matrix  makes  the  plate  bending  out-of-plate  immediately  when  the  plate  is  heated 
without  prebuckling  stage.  The  most  interesting  thing  in  this  study  is  the  mode  shape 
change  and  the  described  solution  procedure  can  automatically  obtain  the  mode  changing 
of  postbuckling  deflection  as  long  as  the  incremental  step  of  temperature  change  is  small 
enough  regardless  of  the  presence  of  mechanical  load.  This  study  demonstrated  that  the 
finite  element  method  can  deal  with  different  planforms,  various  boundary  conditions  as 
expected. 


4.3  Nonlinear  Random  Response 


The  numerical  results  of  random  response  to  acoustic  pressure  only  are  presented  in 
this  section.  A  eight-layer  rectangular  laminate  (  15x  12  in.)  with  the  stacking  sequence 
of  (0/45/-45/90)s  is  considered  first.  The  plate  is  clamped  at  all  four  edges  and  with 
immovable  inplane  boundary  conditions.  Hereafter  it  is  referred  to  as  Panel  1.  The 
material  properties,  mass  density,  and  damping  ratio  are  taken  as: 


Young’s  moduli 
Shear  moduli 
Poisson’s  ratio 
mass  density 
damping  ratio 


El  =  22.5  X  10®,  ^2  =  1.17  X  10®psi 
G23  =  0.4  X  10®,  Gi2  =  Gi3  =  0.66  x  10®psi 
1^12  =  0.22 

p  =  1.45  X  10“'*/6  —  sec? fin.^ 

C  =  0.02 


The  root  mean  square  (RMS)  of  the  maximum  deflection  response  to  a  normal  incidence 
is  shown  in  Fig.  4.13.  The  corresponding  RMS  maximum  micro  strain  is  shown  in  Fig. 
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4.14.  The  response  consists  of  symmetrical  modes,  for  this  laminate  the  2nd,  3rd  and 
4th  modes  are  anti-symmetric,  they  are  not  appeared  in  the  response.  The  frequencies 
are  shown  in  Table  4.3.  Figure  4.15  shows  the  spectrum  density  distribution  of  response 
vs.  frequencies.  When  the  applied  acoustic  load  is  low,  the  response  is  linear.  The  peak 
of  curve  is  very  close  to  the  natural  frequencies.  But  at  high  pressure  levels  the  peaks 
are  shifted  up,  when  the  SPL  reaches  130  dB,  the  first  peak  appears  at  358  Hz. 

The  second  example  is  the  same  as  the  first  one  except  the  stacking  sequence,  it 
is  a  (0/90)  two-layer  laminate.  This  panel  is  referred  to  as  Panel  2.  The  root  mean 
square  of  the  maximum  response  to  a  normal  incidence  is  shown  in  Fig.  4.16.  The 
corresponding  RMS  micro  strain  is  shown  in  Fig.  4.17.  Fig.  4.18  shows  the  spectrum 
density  distribution  of  response  vs.  frequencies. 

The  third  example  is  a  swept  rectangular  plate  as  shown  in  Fig.  4.19.  The  plate 
consists  of  eight-layer  laminate  as  same  as  the  fibrst  one.  It  is  referred  as  Panel  3.  The 
root  mean  square  of  the  maximum  response  to  a  normal  incidence  is  shown  in  Fig.  4.20. 
The  corresponding  RMS  micro  strain  is  shown  in  Fig.  4.21.  Figure  4.22  shows  the 
spectrum  density  distribution  of  response  vs.  frequencies. 

4.4  Acoustic-Thermal  Response 

The  numerical  results  presented  in  this  section  concern  with  the  response  to  combined 
acoustic  and  thermal  loads. 

4.4.1.  Effect  of  Number  of  Modes 

In  order  to  evaluate  the  convergence  characteristics  of  the  present  modal  analysis 
formulation  and  determine  the  required  number  of  modes  for  reasonable  accuracy,  an 
eight-layer  Gr/Ep  rectangular  laminate  (15  x  12  in.)  with  the  stacking  sequence  of  (0/45/- 
45/90)s  is  analyzed  by  mode  number  N=l,  N=2,  N-3  and  N=4.  Anti-symmetrical  modes 
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are  not  included,  The  plate  is  clamped  at  all  four  edges  and  with  immovable  inplane 
boundary  conditions  (u=v=0  at  all  four  edges).  The  full  plate  is  modeled  with  8x8x2 
mesh  or  128  MINS  elements.  For  convenience,  this  plate  is  referred  to  as  the  baseline 
configuration.  The  material  properties,  mass  density,  and  damping  ratio  are  taken  as: 


Young’s  moduli 
Shear  moduli 
Poisson’s  ratio 
mass  density 
damping  ratio 


El  =  22.5  X  10^, E2  =  1.17  X  lO^psi 
G23  =  0.4  X  10^  Gi2  =  Gi3  =  0.66  x  10®psi 
V12  =  0.22 

p  =  1.45  X  10“'^/6  —  se<? 

C  =  0.02 


The  result  of  RMS(Wmax/h)  without  temperature  is  shown  in  Table  4.4.  The  first  12 
mode  characteristics  are  (1,1),  (1,2),  (2,1),  (2,2),  (1,3),  (3,1),  (2,3),  (3,2),  (1,4),  (4,1), 
(4,3)  and  (3,3).  The  first  12  frequencies  are  shown  in  Table  4.5.  The  critical  buckling 
temperature  change  is  37.38®F.  The  result  shows  that  use  of  three  modes  obtained  satisfied 
displacement  results.  Therefore,  in  the  following  calculation,  three  modes  are  used. 
4.4.2.  Effect  of  Thermal  Load 

The  above  configuration  is  analyzed  again  with  ATfATcr  =0.0,  2.0,  3.0.  For 
simplicity,  this  panel  is  called  the  baseline  configuration.  The  root  mean  square  (RMS) 
of  the  maximum  deflection  response  to  a  normal  incidence  and  temperature  change  is 
shown  in  Fig.  4.23.  The  corresponding  RMS  maximum  micro  strain  is  shown  in  Fig. 
4.24.  The  response  consists  of  symmetrical  modes.  For  this  laminate  the  2nd,  3rd  and 
4th  modes  are  antisymmetric,  and  they  do  not  appear  in  the  response.  The  frequencies 
are  shown  in  Table  4.4.  Figure  4.25  shows  the  spectrum  density  distribution  of  response 
vs.  frequencies  at  ATJATcr  =0.0  and  3.0  and  130  dB.  When  the  applied  acoustic  load 
is  low,  the  response  is  linear.  The  thermal  loads  increase  the  nonlinear  stiffness  and  the 
response  is  reduced  due  to  the  large  temperature  rise.  The  peak  of  the  curve  is  very  close 
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to  the  natural  frequencies.  But  at  high  pressure  levels  the  peaks  are  shifted  up,  when  the 
SPL  reaches  130  dB  and  AT/ATcr  =  3.0,  the  first  peak  appears  at  2248  rad/sec. 

For  the  simply  supported  boundary  condition,  the  displacement  response  is  much 
larger  than  clamped  case,  but  the  strain  is  kept  at  about  the  same  amount,  the  results 
for  the  simply  supported  condition  are  shown  in  Figs.  4.26,  4.27  and  4.28.  The  strain 
with  temperature  (Figs.  4.24  and  4.27)  could  be  smaller  or  larger  than  the  one  without 
temperature.  It  illustrates  that  two  effects  are  occurred.  The  thermal  postbuckling 
increases  the  nonlinear  stiffness  which  reduces  the  RMS  deflection,  the  strain  component 
due  to  the  RMS  deflection  is  thus  also  reduced.  On  the  other  hand,  the  thermal  strain 
increases  the  strain  component. 

4,4.3  Antisymmetric  Cross-ply  laminate 

The  material  properties  and  dimension  of  this  clamped  antisymmetric  cross-ply 
laminate  (0/90)  are  the  same  as  the  baseline  configuration.  Figures  4.29  4.30  and  4.31 
show  the  mode  shapes  of  this  panel  at  AT=0, 97.972  and  149.688°F,  respectively.  In  these 
figures  two  features  should  be  noticed:  the  mode  sequence  is  changed  with  temperature 
rise;  and  some  mode  shapes  are  not  exactly  symmetric  or  antisymmetric  as  in  the  case 
of  isotropic  material.  Figure  4.32  shows  the  displacement  response  while  Figure  4.33 
shows  the  micro  strain  distribution  vs.  sound  pressure  level. 
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4.4.4.  Skewed  Panel  and  Non-uniform  Temperature 


The  planform  of  this  panel  is  shown  in  Fig.  4.19.  the  material  property  is  the  same 
as  the  baseline  configuration.  The  non-uniform  temperature  change  is  that  at  the  edge 
grid  points  the  temperature  change  is  zero,  at  interior  grids  is  uniform.  The  boundary 
condition  is  clamped.  The  results  are  shown  in  Figs.  4.34  and  4.35.  The  critical  buckling 
temperature  changes  for  uniform  and  non-uniform  temperature  distributions  are  67.36°F 
and  76.32°F,  respectively.  It  can  be  seen  from  the  figures  that  the  responses  and  strains 
for  uniform  temperature  and  nonuniform  temperature  have  very  little  difference.  This 
illustrates  that  the  temperature  gradient  along  the  edge  has  little  influence  on  random 
responses. 

4.4.5.  Grazing  Incidence  Wave 

The  plate  studied  in  this  case  is  the  same  as  the  baseline  configuration.  But  the 
boundary  condition  is  simply  supported  for  transverse  displacement  and  immovable  for 
in-plane  displacements.  The  thickness  of  the  plate  is  0.048  in.  The  dimension  is  15  in. 
by  12  in.  The  stacking  sequence  is  (0/45/-45/90)s. 

The  result  of  this  example  is  very  interesting.  Because  the  acoustic  wave  is  travelling 
along  the  positive  direction  of  x  axis  with  a  speed  c,  the  acoustic  pressure  on  the  plate 
along  the  x-axis  is  not  uniform  and  the  antisymmetric  modes  about  y-axis  participates 
in  the  response  of  the  plate.  Therefore  the  maximum  deflection  point  moves  forward 
slightly  as  shown  in  Fig.  4.36.  The  deflection  spectrum  density  is  shown  in  Fig.  4.37. 

4.4.6  Effect  of  Initial  Imperfection 

If  the  plate  has  some  initial  imperfection  in  deflection,  the  nonlinear  stiffness  due  to 
initial  deflection  reduces  the  response  as  compared  to  flat  plate.  It  is  also  stiffer  than 
thermal  postbuckling  deflection,  assuming  that  they  have  the  same  maximum  deflection 
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as  shown  in  Fig.  4.26.  For  thermal  postbuckling,  the  panel  is  also  subjected  to  certain 
thermal  stress.  For  initial  deflection,  the  plate  has  only  geometric  stiffness  which  reduces 
the  random  response;  while  thermal  postbuckling  plate  has  thermal  stresses  and  thermal 
deflection.  The  results  are  shown  in  Figs.  4.38,  4.39  and  4.40. 
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Table  4.1  Results  for  a  rhombic  cantilever  plate 


Deflection  at  locations  (in.) 

1 

2 

3 

4 

5 

6 

MIN  3  result 

0.263 

0.178 

0.108 

0.103 

0.048 

0.019 

Experiment  value 

0.297 

0.204 

0.121 

0.129 

0.056 

0.022 

Table  4.2  Results  of  RMS  (Wmax/h)  and  Micro-strain  for  simply  supported  plate 

(N=number  of  modes) 


Present:N=l 

resent:N=2 

present:  N=4 

Chiang:N=l 

Chiang:N=4 

SSL  ‘ 

RMS 

Micro 

RMS 

Micro 

RMS 

Micro 

RMS 

Micro 

RMS 

Micro 

(dB) 

W 

strain 

W 

strain 

W 

strain 

W 

strain 

W 

strain 

ig 

1.04 

64.5 

1.04 

80.5 

1.07 

77.7 

1.030 

87.4 

1.031 

112.0 

H 

1.89 

213.2 

1.89 

274.5 

1.96 

259.2 

1.902 

256.5 

1.905 

361.0 

Table  4.3  Linear  Natural  Frequencies(Hz)  of  Panel  1 


1st 

2nd 

3rd 

4th 

5th 

6th 

7th 

100.3 

184.9 

210.4 

286.0 

321.0 

360.8 

413.6 

The  2nd,  3rd  and  4th  are  antisymmetrical  modes. 
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Table  4.4  Convergence  of  RMS(Wmax/h)  with  Number  of  the  Modes 


No.  of  Modes 

N=1 

N=2 

N=3 

N=4 

RMS(Wmax/h) 

2.0252 

1.9626 

1.8676 

1.8695 

Table  4.5  The  Natural  Frequencies  (rad/sec.) 


1st 

2nd 

3rd 

4th 

5th 

6th 

629.95 

1161.6 

1796.8 

7th 

8th 

9th 

10th 

11th 

12th 

2598.7 

2758.1 

3125.7 

3393.4 

3526.6 

3734.8 

Table  4.6  Frequencies  for  tiie  baseline  configuration  (rad/sec.) 


AT/ATcr 

1st 

2nd 

3rd 

4th 

5th 

6th 

— 

0.0 

630.0 

1161.6 

1322.0 

1796.8 

2016.6 

2267.0 

2.0 

811.9 

848.3 

1135.6 

1406.0 

1784.2 

2122.8 

3.0 

1046.2 

1100.7 

1421.6 

1940.2 

2251.4 
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Fig.  4.1  Twisted  ribbon  test  via  transverse  comer  forces 
(E=10^,  i/=0.25,  h=0.05)(From  Tessler  and  Hughes,  1985) 


Fig.  4.2  Twisted  ribbon  test  via  twisting  comer  moments(E=10^, 
i^=0.25,  h=0.05)(From  Tessler  and  Hughes,  1985) 
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Wmax/h 


q=0.26066  psi 
E=10.5  E06  psi 
1=0.125  in. 
Nu=0.3 

y  4x4x2  mesh 


Fig.  4.3  Rhombic  cantilever  plate 


Fig.  4.4  Wmax/h  vs.  AT/ATref.  for  (0/90)  composite  plate 
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AT(deg.  of  F) 


Fig.  4.5  Wmax/h  vs.  AT  for  an  isosceles  triangular  plate 


0  100  200  300  400 


Length  to  thickness  ratio  (a/h) 


Fig.  4.6  Dimensionless  critical  temperature  vs.  a/h  for  a  rectangular  laminate 
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''0  80  90  100 

AT 

Fig.  4.7  Wmax/h  vs.  AT  for  various  number  of  layers 
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Fig.  4.9  Centerline  deflection  of  a  36x12x0.048  in.  angle-ply  (60/-60)  laminate 


Fig.  4.10  Wmax/h  vs.  AT  for  a  long  rectangular 
laminate  with  mechanical  load  0.01  psi 
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W/h 


Fig.  4.11  Centerline  deflection  of  a  36x12x0.048  in. 
angle-ply  (60/-60)  laminate  with  mechanical  load  0.01  psi 


Fig.  4.12  Wmax/h  vs.  skew  angle  for  rectangular  plates 
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RMS  max.  micro  strain  RMS  Wmax/h 


Sound  Spectrum  Level  (dB) 

Fig.  4.13  RMS  Wmax/h  vs.  SPL  for  Panel  1 


Fig.  4.14  RMS  max.  micro  strain  vs.  SPL  for  Panel  1 
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Fig.  4.15  The  maximum  deflection  spectmm  vs.  frequency  for  Panel  1 


(0/90),  15”X12"X0.048" 
Qamped 


Deflection  Spectrum  Gw 


Fig.  4.18  The  maximum  deflection  spectmm  vs.  frequency  for  Panel  2 
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Fig.  4.19  Planform  of  Panel  3 


Fig.  4.20  RMS  Wmax/h  vs.  SPL  for  Panel  3 
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Fig.  4.21  RMS  max.  micro  strain  vs.  SPL  for  Panel  3 


Fig.  4.22  The  maximum  deflection  spectrum  vs.  frequency  for  Panel  3 
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Fig.  4.23  RMS  (Wmax/h)  vs.  SPL  for  the  baseline  configuration 


Fig.  4.24  Micro-strain  vs.  SPL  for  the  baseline  configuration 
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Deflection  Spectrum  Gw 


Fig.  4.26  RMS  (Wmax/h)  vs.  SPL  for  the  simply  supported  panel 


Daflectfon  Spectrum  Gw 


AT/ATcr=0 

AT/ATcrs2 

AT/ATcr^ 


(0/45M5/0)s 
Simply  supported 
15X12X0.048  in. 
128  element.  N=3 


0  0 


Fig.  4.29  The  mode  shapes  of  (0/90)  clamped  panel 


2nd  mode 


0  0 

Fig.  4.29  Continued 
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3rd  mode 


0  0 


Fig.  4.29  Continued 


4th  mode 


0  0 


Fig.  4.29  Continued 


0  0 


Fig.  4.29  Continued 

6th  mode 
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0  0 

Fig.  4.30  The  mode  shapes  of  (0/90)  clamped  panel 
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2nd  mode,  Tern.  Changes97.792  (deg.  F) 


Fig.  4.30  Continued 


3rd  mode,  Tern.  Change=97.792  (deg.  F) 


0  0 


Fig.  4.30  Continued 


4th  mode,  Tern.  Chang©=:97.792  (deg.  F) 


0  0 


Fig.  4.30  Continued 


5th  mode.  Tern.  Change=97.792  (deg.  F) 


Fig.  4.30  Continued 
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6th  mode,  Tern.  Change=97.792  {deg.  F) 


0  0 

Fig.  4.30  Continued 


Fig.  4.31  The  mode  shapes  of  (0/90)  clamped  panel 


2nd  mode,  Tern.  Change^l  49.688  (deg.  F) 


0  0 

Fig.  4.31  Continued 


3rd  mode,  Tem.  Change* 149.688  (deg.  F) 


0  0 

Fig.  4.31  Continued 
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5lh  mode,  Tern.  Change^  149.688  (deg.  F) 


0  0 

Fig.  4.31  Continued 
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1000 


Fig.  4.33  RMS  max.  micro  strain  vs.  SPL  for  the  (0/90)  panel 


Fig.  4.34  RMS  (Wmax/h)  vs.  SPL  for  a  skewed  panel 


’9 


Fig.  4.35  RMS  max.  micro  strain  vs.  SPL  for  a  skewed  panel 


Fig.  4.36  Distributions  of  RMS  W/h  along  the  center  line  of  the  panel 
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RMS  (Wmax/h)  Deflection  Spectrum  Gw 


Fig.  4.37  The  maximum  deflection  spectram  vs.  frequency 
for  the  simply  supported  panel  subjected  to  grazing  wave 


Sound  Pressure  Level  (dB) 


Fig.  4.38  RMS  (Wmax/h)  vs.  SPL  for  the 
simply  supported  panel  with  initial  deflection  Wo 


RMS  max.  miao-strain  RMS  max.  micro-strain 


Fig.  4.39  RMS  max.  Micro-strain  vs.  SPL  for  the 
simply  supported  plate  with  initial  deflection  Wo 


Fig.  4.40  RMS  max.  Micro-strain  vs.  SPL  for  the 
simply  supported  plate  with  initial  imperfection  Wo 
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Chapter  5 


CONCLUSIONS 

Using  three-node  Mindlin  plate  element  with  improved  shear  deformation,  the  gov¬ 
erning  nonlinear  equations  of  motion  have  been  derived  for  composite  structures  subjected 
to  a  combined  acoustic/thermal  loading.  In  order  to  simulate  the  acoustic  waves  in  the 
progressive  wave  test  facility,  a  grazing  incidence  wave  model  is  used  in  the  derivation. 
An  innovative  solution  procedure  has  been  created  and  the  equations  of  motion  were 
solved  for  applications  of  thermal  postbucking  and  large  deflection  random  response  of 
thermally  buckled  structures. 

The  critical  temperature  change  that  produced  panel  buckling  was  determined  by  the 
incremental  equations  of  motion.  The  first  buckling  mode  was  used  as  the  initial  shape 
of  the  postbuckling  solution.  Newton-Raphson  iteration  method  was  used  to  solve  for  the 
deflections  corresponding  to  a  given  temperature  raise  distribution.  The  extension  and 
bending  coupling  makes  the  plate  bending  out  of  plane  immediately  when  the  plate  is 
heated  without  prebuckling  stage.  The  most  interesting  aspect  of  this  study  is  the  mode 
shape  change.  The  described  solution  procedure  can  automatically  obtains  the  mode 
changing  in  the  postbuckling  stage  as  long  as  the  incremental  step  of  temperature  change 
is  small  enough,  regardless  of  the  presence  of  mechanical  load. 

In  order  to  solve  the  system  equation  of  motion,  an  innovative  solution  procedure 
is  described.  The  response  is  assumed  to  be  the  sum  of  static  component  and  dynamic 
component.  Substituting  the  total  displacements  into  the  system  equation  of  motion 
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and  regrouping  the  terms,  the  equation  is  of  the  form  =  G({W}J.  The 

left-hand-side  of  the  equation,  ■^({VF}^),  is  independent  of  time  t;  while  the  right-hand- 
side,  (j({W}J,  is  time  dependent.  Therefore,  the  only  possibility  for  both  F  and  G  to 
exist  is  that  both  F  and  G  equal  to  zero,  and  the  two  equations  are  thus  obtained.  For 
the  dynamic  components,  the  modal  transformation  was  used  to  reduce  the  number  of 
equations.  Then  the  equivalent  linearization  method  is  utilized  to  obtain  the  coupled  linear 
equation.  Finally,  the  modal  transform  was  used  to  uncouple  the  equation  of  motion,  and 
an  iteration  was  used  to  obtain  the  random  response. 

The  most  significant  contributions  of  the  present  study  are  the  formulation  and 
solution  procedure,  including  grazing  wave,  of  nonlinear  modal  equations  used  to  describe 
the  random  response  of  composite  structures  to  combined  acoustic  and  thermal  loads. 
These  general  equations  are  applicable  not  only  to  the  present  research,  but  also  to  other 
dynamic  problems  like  gust  response,  buffet  response  of  an  aircraft.  Comparing  with  the 
Sequential  Load  method(SQL),  this  SiMultaneous  I^ad  method(SML)  is  much  logical 
mathematically  and  straightforward,  it  is  then  easier  to  formulate  nonlinear  problems  with 
combined  loading.  The  solution  procedure  itself  can  take  care  of  the  inter-dependence 
between  the  thermal  effects  and  the  acoustic-structural  response.  In  the  SQL  method 
engineering  judgment  is  essential,  otherwise  some  terms  might  be  missed. 

From  the  results,  an  interesting  observation  is  that  the  antisymmetric  modes  partici¬ 
pate  in  the  response  of  the  plate  for  grazing  incidence  acoustic  wave.  It  is  demonstrated 
that  three  or  four  modes  will  give  converged  RMS  deflection.  It  is  also  found  that 
the  RMS  maximum  strain  with  temperature  could  be  either  smaller  or  larger  than  the 
one  without  temperature.  This  is  due  to  that:  (1)  the  temperature  increases  the  thermal 
component,  and  (2)  the  thermal  postbuckling  deflection  increases  the  nonlinear  stifftiess 
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which  reduces  the  RMS  deflection  and  it  leads  to  smaller  strain  component.  Uniform 
and  nonuniform  temperature  distribution  effect  on  random  responses  is  investigated.  The 
nonuniform  temperature  considered  is  that  there  is  a  temperature  gradient  along  the  edge 
of  the  plate.  The  results  show  that  there  is  very  little  difference  in  random  responses  for 
the  two  temperature  distributions  studied.  For  plate  with  initial  imperfection  in  deflection 
which  has  the  same  maximum  deflection  as  the  thermal  postbuckling  deflection,  the  plate 
with  initial  imperfection  is  stififer  and  leads  to  smaller  random  responses. 

Future  improvements  for  grazing  wave  model  and  numerical  integration  methods  are 
needed.  The  correlation  study  between  numerical  results  and  test  results  is  very  important 
if  one  wants  to  use  this  method  to  a  practical  problem. 
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Appendix  A 

THE  FORMULATIONS  FOR 


Cij,fn  =  2  (^^126  -  2/3i6)y23 


(A.1) 


=  ^(?/236  -  2/126)2/31 


(A.2) 


C'V’^13  =  2  (^^316  ~  2/236)2/12 


(A.3) 


Crj,^u  =  ^(2:136  -  ®2i6)2/23  -  \x2ii\yzi  +  -^xiziiyvi 


(A.4) 


C'V’^/’lS  =  ^(®2i6  -  2:326)2/31  -  ^2:3262/12  +  2^2162/23 


(A.5) 


C^^ie  =  ^(2:326  -  2:136)2/12  -  ^2:1362/23  +  2^3262/31 


(A.6) 


C'v,^2i  =  ^(2/126  -  2/316)2:32  —  2^2162/31  +  2^1362/12  (A.7) 


1  11 

^*^^22  =  2(^/236  -  2/126)2:13  —  2^3262/12  +  2^2162/23 


(A.8) 


1  11 

C'^^/»23  =  2(2/316  -  2/236)2:21  -  2^1362/23  +  2^3262/31 


(A.9) 
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(A.10) 


=  2  (^136  “  2:216)2:32 


C'^^&25  =  2  (^216  -  2:326)2:13 


C'^/>^26  =  2(^326  -  2:136)2:21 


(A.11) 


(A.12) 
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Appendix  B 

THE  ELEMENT  MATRICES 


Linear  stiffness  matrix 

W*  =  j  [Cif[D]mdA  (B.l) 

A 

=  J  nf[B][C„]dA  (B.2) 

A 

=  /  [Cr^flBmdA  (B.3) 

A 

[k]r.=  j\Cmf[A]\Cm\dA  (B.4) 

A 

Linear  stiffness  matrix  due  to  Wo(x,y) 

[Mi=  /  \Ctt]^\<lof[Am[Ci,t]dA  (B.5) 

A 

[k.]H  =  /  [CiA\hf[B\mdAA  I  [C^,{f[e/\Am[C^,^]dA  (B.6) 

^  A 

Wto  =  /  [Ci,if[dof\A][C„]dA  (B.7) 

A 

=  /  [Clf[Bm[Ci,i]dAA  I  [C^^]'^[»„f[A]m[C^i]dA  (B.8) 

A  A 

fe)*=  f  lCti.f[eofmCt]dA+  f  \Cif\B]%\[Ci,i,]dA 

'I  (B.9) 

+  /  [Ci^]  \Kf[A][h\\Ci,t\dA 

A 

%\^  =  /  [CftY[h?[A][Cr,idA  (B.IO) 

A 
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—  J  [Cm]  [-^][^o]  [C^ilf^ldA  (B.ll) 

A 

Linear  stiffness  matrix  due  to  {N^^} 

[kNAT]h  =  J  [C^h]'^[NAT][C^j;h]dA  (B.12) 

[kNAT]h^  =  j  (B.13) 

A 

[kNAT]^h  =  J  [CrPi;]'^[NAT][C^h]dA  (B.14) 

A 

[%at]^  =  j  (B.15) 

A 

Linear  stiffness  matrices  due  to  {Nq} 

[kNo]i^  J  [C^bf[No][C^i]dA  (B.16) 

A 

[kNo]b^  =  j  [C^hf[No][C^^]dA  (B.17) 

A 

[^iVo]^j  “  J  [-^o]  (B.18) 

A 

[kNo]^  =  J  [C^i;f[No][C^^]dA  (B.19) 

A 

First-order  nonlinear  stiffness  matrix 

I"llw  =  /  [Ci,if[efmCi]dA  (B.20) 

A 

["IJim  =  J  [C^Af[ef[A][C^]dA  (B.21) 

A 

[nlj^j  =  j  [Cif[B]\e][C^i]dA  (B.22) 

A 

["!)*  =  J  [C^f[Bm[C^,^,]dAA  f  [C^^Y[6f[B][Ct]dA  (B.23) 

A  A 
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(B.24) 


['>!)♦„  =  /  [C^^YW[AlC„\dA 

A 

=  /  {Cj^[Am[C^,i]iA  (B.25) 

A 

[nl]„^  =  j  [C^f[Am  [Ct^]dA  (B.26) 

A 

First-order  nonlinear  stiffiiess  matrix  due  to  Wo(x,y) 

[nl.],  =  J  [C^if[ef[A][e,][C^4dA  +  j  [C^t]^[«,f[Am[C^,i]dA  (B.27) 

A  A 

[«1<.W  =  /  [C^hf[ef[A\%][Cn,]dA  +  j  [C^,if[6.f[A]m[C^,^,]dA  (B.28) 

A  A 

[«1.U=  /  [Ci,*f\»f\Am[C^,i]dA  +  I  [C^,^Y[e,f[A]\()][C^,^]dA  (B.29) 

A  A 

[nlo\^  =  J  [C^^f[ef[A][eo][C^^]dA  +  J  [C^^f[9of[A][0][C^^]dA  (B.30) 

A  A 

First-order  nonlinear  stififhess  matrix  due  to  {Nni}(=[A]{e^}) 

[nljVmlft  =  I  [C,l;bf[Nm][C^b]dA  (B.31) 

A 

[nliv^ls#  =  j  [c^s]  ^IiV„l  [C**]  dA  (B.32) 

A 

[nlNmU  =  J  [Ci,^f[Nm][C^b]dA  (B.33) 

A 

=  y*  [C,i;4,]^[Nm\[Cji,^]dA  (B.34) 

A 

First-order  nonlinear  stiffness  matrix  due  to  {NbK=[B]{K}) 

[i^'^Nblb  =  J  [^ipb]^Wb][0^b]dA  (B.35) 

A 

[f^'^Nblb^  =  J  [C^b]^Wb][C^ik\^A  (B.36) 

A 
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MjVi]v.6  =  J  [Cf^]^[Ni][C^i]dA 

A 

(B.37) 

[nlNb]^=  J  lC^^fm[C^^]dA 

A 

(B.38) 

Second-order  nonlinear  stifi&iess  matrix 

[n21,  =  11  [C^^]^[6f[Am[C^,i\dA 

A 

(B.39) 

["21w  =  1/  [C^,^flef^AM[C^,^,]dA 

A 

(B.40) 

["21#  =  1/  [C^,tf[ef[Am[Cti]dA 

A 

(B.41) 

1"2]^  =  [C^,^f\»f\A]\e][C^^]dA 

A 

(B.42) 

Linear  stiffness  matrix  due  to  shear 

[^.lt  =  j  [C^i\^[A.\[C,t]dA 

A 

(B.43) 

[^5]^^  =  J  [Cyh]  [^s] 

A 

(B.44) 

[^s]^h  ~  J  [^7^]  [^3]  [Cjb]  dA 

A 

(B.45) 

l^s]^  =  J  [C^7^]  [-^s] 

A 

(B.46) 

Load  vectors 

{PAT}m  =  J  [Cmf  {N^T}dA 

A 

(B.47) 

{PNo}m  =  J -[Cmf{No}dA 

A 

(B.48) 

{PATo}b  =  J  [0^b]'^[^of' {N/:^T}dA 

A 

(B.49) 
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{PN>oh  =  -f  [C^,if{Dof{Nc]dA 

A 

{pATo}*  =  J  [C^4fVof{NAT]dA 

A 

{pno.}^,  =  -[  [C^^f[hf{No}dA 

A 

{pat}*  =  J  [Cif  {Mcjr)dA 

A 

{Pp)h  =  J  [Hwf  p{xyy,t)dA 

A 

OO 

=  ~  J  j  dAdu 

“OO  A 

OO 

=  i  /  P{u>m^)}e-“‘d^ 

— OO 


=  ^  J  P{iv)e^*  J  e-^^^^[H^^fdAdu; 

-OO  A 


OO 

—  OO 


(  For  Eqs.  (2.107)  and  (2.108)  see  Ref. [Clarkson]) 


Mass  matrix 


[m\  =  J  [Hy,f  ph[H^]dA 

A 


Nj*  =  J  IH^f  ph[H^^,]dA 

A 


A 


(B.50) 

(B.51) 

(B.52) 

(B.53) 

(B.54) 


(B.55) 


(B.56) 

(B.57) 

(B.58) 
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LIST  OF  SYMBOLS 


a,  b 

[A],[B],[D] 

c 

[C] 

E 

E[] 

{e} 

{/},  {/} 

G 

Hjiu) 

[H\ 

h 

[K] 

{M} 

[M] 

{N} 

[Nl]  ,[N2] 

{P} 

{q} 

[Q] 

{R} 

t 

u,  V 
w 

x,y,z 

Greeek  Symbols 
a 

AT 

uy,  {7} 


lengths  of  a  plate 

laminate  extensional,  extension-bending  and  bending  stiffnesses 
acoustic  wave  traveling  speed 

the  matrices  relating  slope  and  curvature  with  displacements 
Young’s  modulus 
expected  value 
linearization  error  vector 

modal  force  vector  of  coupled  and  uncoupled  modal  equations 

the  stiffness  term  of  the  modal  equations  of  motion 

shear  modulus  of  the  materials 

single  degree  of  freedom  transfer  function 

shape  function  matrices 

thickness  of  plate 

linear  stiffness 

resultant  moments  per  unit  length  acting  on  a  laminate 
mass  matrix 

resultant  forces  per  unit  lengtii  acting  on  a  laminate 

die  first  and  second  order  nonlinear  system  stiffness  matrices 

nodal  load  vector 

coupled  modal  coordinate  vector 

lamina  reduced  stiffness 

resultant  shear  forces  per  unit  lengfli  acting  on  a  laminate 
time 

inplane  displacement 
lateral  displacement 
Cartesian  coordinates 

thermal  expansion  coefficient  and  shear  correction  factor 
temperature  distribution 
strain  vector 
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c 

damping  ratio 

n 

uncoupled  modal  coordinates 

em 

acoustic  wave  incident  angle,  and  as  defined  in  Eq.  (2.26) 

[<!>] 

modal  matrix 

fpxyfky 

rotations  of  the  normal  around  x  and  y  axes  due  to  bending  only 

Q 

equivalent  linear  frequency 

Subscripts 

linear  frequency  of  the  deformed  structure 

b 

bending 

m 

membrane 

Nb 

stiffness  matrices  due  to  {Nb} 

Nm 

stiffness  matrices  due  to  {Nm} 

No 

stiffness  matrices  due  to  {N©} 

NAT 

stiffness  matrices  due  to  {N^^} 

0 

quantity  related  with  initial  displacement 

T 

tangent  stiffriess 

AT 

thermal 
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